H₂SO₄–H₂O–NH₃ ternary ion-mediated nucleation (TIMN): kinetic-based model and comparison with CLOUD measurements

New particle formation (NPF) is known to be an important source of atmospheric particles that impacts air quality, hydrological cycle, and climate. Although laboratory measurements indicate that ammonia enhances NPF, the physicochemical processes underlying the observed effect of ammonia on NPF are yet to be understood. Here we present a comprehensive kinetically based H₂SO₄–H₂O–NH₃ ternary ion-mediated nucleation (TIMN) model that is based on the thermodynamic data derived from both quantum-chemical calculations and laboratory measurements. NH₃ was found to reduce nucleation barriers for neutral, positively charged, and negatively charged clusters differently, due to large differences in the binding strength of NH₃, H₂O, and H₂SO₄ to small clusters of different charging states. The model reveals the general favor of nucleation of negative ions, followed by nucleation on positive ions and neutral nucleation, for which higher NH₃ concentrations are needed, in excellent agreement with Cosmics Leaving OUtdoor Droplets (CLOUD) measurements. The TIMN model explicitly resolves dependences of nucleation rates on all the key controlling parameters and captures the absolute values of nucleation rates as well as the dependence of TIMN rates on concentrations of NH₃ and H₂SO₄, ionization rates, temperature, and relative humidity observed in the well-controlled CLOUD measurements well. The kinetic model offers physicochemical insights into the ternary nucleation process and provides a physics-based approach to calculate TIMN rates under a wide range of atmospheric conditions.

1 Introduction

New particle formation (NPF), an important source of particles in the atmosphere, is a dynamic process involving interactions among precursor gas molecules, small clusters, and preexisting particles (Yu and Turco, 2001; Zhang et al., 2012). H₂SO₄ and H₂O are known to play an important role in atmospheric particle formation (e.g., Doyle, 1961). In typical atmospheric conditions, the species dominating the formation and growth of small clusters is H₂SO₄. The contribution of H₂O to the nucleation is related to the hydration of H₂SO₄ clusters (or, in the other words, modification of the composition of nucleating clusters), which reduces the H₂SO₄ vapor pressure and hence diminishes the evaporation of H₂SO₄ from the pre-nucleation clusters. NH₃, the most abundant gas-phase base molecule in the atmosphere and a very efficient neutralizer of sulfuric acid solutions, has long been proposed to enhance nucleation in the lower troposphere (Coffman and Hegg, 1995), although it has been well recognized that earlier versions of the classical ternary nucleation model (Coffman and Hegg, 1995; Korhonen et al., 1999; Napari et al., 2002) significantly overpredict the effect of ammonia (Yu, 2006a; Merikanto et al., 2007; Zhang et al., 2010).

The impacts of NH₃ on NPF have been investigated in a number of laboratory studies (Kim et al., 1998; Ball et al., 1999; Hanson and Eisele, 2002; Benson et al., 2009; Kirkby et al., 2011; Zollner et al., 2012; Froyd and Lovejoy, 2012; Glasoe et al., 2015; Schobesberger et al., 2015; Kürten et al., 2016) including those recently conducted at the European Organization for Nuclear Research (CERN) in the framework of the CLOUD (Cosmics Leaving OUtdoor Droplets) experiment that has provided a unique dataset for quantitatively examining the dependences of ternary H₂SO₄–H₂O–NH₃ nucleation rates on concentrations of NH₃ ([NH₃]) and H₂SO₄ ([H₂SO₄]), ionization rate (Q), temperature (T), and relative humidity (RH) (Kirkby et al., 2011; Kürten et al., 2016). The experimental conditions in the CLOUD chamber, a 26.1 m³ stainless steel cylinder, were well controlled, while impacts of potential contaminants were minimized (Schnitzhofer et al., 2014; Duplissy et al., 2016). Based on CLOUD measurements in H₂SO₄–H₂O–NH₃ vapor mixtures, Kirkby et al. (2011) reported that an increase in [NH₃] from ∼0.03 ppb (parts per billion, by volume) to ∼0.2 ppb can enhance ion-mediated (or induced) nucleation (IMN) rate by 2–3 orders of magnitude and that the IMN rate is a factor of 2 to > 10 higher than that of neutral nucleation under a typical level of contamination by amines. In the presence of ionization, common highly polar atmospheric nucleation precursors such as H₂SO₄, H₂O, and NH₃ molecules tend to cluster around ions, and charged clusters are generally much more stable than their neutral counterparts with enhanced growth rates as a result of dipole–charge interactions (Yu and Turco, 2001).

Despite various laboratory measurements indicating that ammonia enhances NPF, the physicochemical processes underlying the observed different effects of ammonia on the formation of neutral, positively charged, and negatively charged clusters (Schobesberger et al., 2015) are yet to be understood. To achieve such an understanding, a nucleation model based on the first principles is needed. Such a model is also necessary to extrapolate data obtained in a limited number of experimental conditions to a wide range of atmospheric conditions, in which [NH₃], [H₂SO₄], ionization rates, T, RH, and surface areas of preexisting particles vary widely depending on the region, pollution level, and season. The present work aims to address these issues by developing a kinetically based H₂SO₄–H₂O–NH₃ ternary IMN (TIMN) model that is based on the molecular clustering thermodynamic data. The model predictions are compared with relevant CLOUD measurements and previous studies.

2 Kinetic-based H₂SO₄–H₂O–NH₃ ternary ion-mediated nucleation (TIMN) model

2.1 Background

Most nucleation models developed in the past for H₂SO₄–H₂O binary homogeneous nucleation (e.g., Vehkamäki et al., 2002), H₂SO₄–H₂O ion-induced nucleation (IIN; e.g., Hamill et al., 1982; Raes et al., 1986; Laakso et al., 2003), and H₂SO₄–H₂O–NH₃ ternary homogeneous nucleation (Coffman and Hegg, 1995; Korhonen et al., 1999; Napari et al., 2002) have been based on the classical approach, which employs capillarity approximation (i.e., assuming that small clusters have the same properties as bulk) and calculates nucleation rates according to the free energy change associated with the formation of a “critical embryo”. Yu and Turco (1997, 2000, 2001) developed a neutral and charged binary H₂SO₄–H₂O nucleation model using a kinetic approach that explicitly treats the complex interactions among small air ions, neutral and charged clusters of various sizes, precursor vapor molecules, and preexisting aerosols. The formation and evolution of cluster size distributions for positively and negatively charged cluster ions and neutral clusters affected by ionization, recombination, neutralization, condensation, evaporation, coagulation, and scavenging have been named IMN (Yu and Turco, 2000). The IMN theory significantly differs from classical IIN theory (e.g., Hamill et al., 1982; Raes et al., 1986; Laakso et al., 2003), which is based on a simple modification of the free energy for the formation of a critical embryo by including the electrostatic potential energy induced by the embedded charge (i.e., Thomson effect; Thomson, 1888). The classical approach does not properly account for the kinetic limitation to embryo development, enhanced stability and growth of charged clusters associated with dipole–charge interaction (Nadykto and Yu, 2003; Yu, 2005), and the important contribution of neutral clusters resulting from ion–ion recombination to nucleation (Yu and Turco, 2011). In contrast, these important physical processes are explicitly considered in the kinetic-based IMN model (Yu, 2006b).

Since the beginning of the century, nucleation models based on the kinetic approach have also been developed in a number of research groups (Lovejoy et al., 2004; Sorokin et al., 2006; Chen et al., 2012; Dawson et al., 2012; McGrath et al., 2012). Lovejoy et al. (2004) developed a kinetic ion nucleation model, which explicitly treats the evaporation of small neutral and negatively charged H₂SO₄–H₂O clusters. The thermodynamic data used in their model were obtained from measurements of small ion clusters, ab initio calculations, the thermodynamic cycle, and some approximations (adjustment of Gibbs free energy for neutral clusters calculated based on liquid droplet model, interpolation, etc.). Lovejoy et al. (2004) did not consider the nucleation on positive ions. Sorokin et al. (2006) developed an ion cluster–aerosol kinetic (ICAK) model, which uses the thermodynamic data reported in Froyd and Lovejoy (2003a, b) and empirical correction terms proposed by Lovejoy et al. (2004). Sorokin et al. (2006) used the ICAK model to simulate dynamics of neutral and charged H₂SO₄–H₂O cluster formation and compared the modeling results with their laboratory measurements. Chen et al. (2012) developed an approach for modeling NPF based on a sequence of acid–base reactions, with sulfuric acid evaporation rates (from clusters) estimated empirically based on measurements of neutral molecular clusters taken in Mexico City and Atlanta. Dawson et al. (2012) presented a semiempirical kinetics model for nucleation of methanesulfonic acid (MSA), amines, and water that explicitly accounted for the sequence of reactions leading to formation of stable particles. The kinetic models of Chen et al. (2012) and Dawson et al. (2012) consider only neutral clusters.

McGrath et al. (2012) developed the Atmospheric Cluster Dynamics Code (ACDC) to model the cluster kinetics by solving the birth–death equations explicitly, with evaporation rate coefficients derived from formation free energies calculated by quantum chemical methods (Almeida et al., 2013; Olenius et al., 2013). The ACDC model applied to the H₂SO₄–dimethylamine (DMA) system considers zero to four base molecules and zero to four sulfuric acid molecules (Almeida et al., 2013). Olenius et al. (2013) applied the ACDC model to simulate the steady-state concentrations and kinetics of neutral and negatively and positively charged clusters containing up to five H₂SO₄ and five NH₃ molecules. In ACDC, the nucleation rate is calculated as the rate of clusters growing larger than the upper bounds of the simulated system (i.e., clusters containing four or five H₂SO₄ molecules) (Kürten et al., 2016).

The kinetic IMN model developed by Yu and Turco (1997, 2001) explicitly simulates the dynamics of neutral, positively charged, and negatively charged clusters, based on a discrete-sectional bin structure that covers the clusters containing 0, 1, 2,…, 15,…H₂SO₄ molecules to particles containing thousands of H₂SO₄ (and H₂O) molecules. In the first version of the kinetic IMN model (Yu and Turco, 1997, 2001), due to the lack of thermodynamic data for the small clusters, the compositions of neutral and charged clusters were assumed to be the same and the evaporation of small clusters was accounted for using a simple adjustment to the condensation accommodation coefficients. Yu (2006b) developed a second-generation IMN model that incorporated newer thermodynamic data (Froyd, 2002; Wilhelm et al., 2004) and physical algorithms (Froyd, 2002; Wilhelm et al., 2004) and explicitly treated the evaporation of neutral and charged clusters. Yu (2007) further improved the IMN model by using two independent measurements (Marti et al., 1997; Hanson and Eisele, 2000) to constrain monomer hydration in the H₂SO₄–H₂O system and by incorporating experimentally determined energetics of small neutral H₂SO₄–H₂O clusters that became available then (Hanson and Lovejoy, 2006; Kazil et al., 2007). The first and second generations of the IMN model were developed for the H₂SO₄–H₂O binary system, although the possible effects of ternary species such as the impact of NH₃ on the stability of both neutral and charged pre-nucleation clusters have been pointed out in these previous studies (Yu and Turco, 2001; Yu, 2006b). The present work extends the previous versions of the IMN model in the binary H₂SO₄–H₂O system to the ternary H₂SO₄–H₂O–NH₃ system, as described below.

Figure 1 Schematic illustration of kinetic processes controlling the evolution of positively charged (H⁺S_sW_wA_a), neutral (S_sW_wA_a), and negatively charged (${\mathrm{S}}^{-}{\mathrm{S}}_{s-\mathrm{1}}{\mathrm{W}}_{\mathrm{w}}{\mathrm{A}}_{\mathrm{a}}$) clusters–droplets that are explicitly simulated in the ternary ion-mediated nucleation (TIMN) model. Here S, W, and A represent sulfuric acid (H₂SO₄), water (H₂O), and ammonia (NH₃), respectively, while s, w, and a refer to the number of S, W, and A molecules in the clusters–droplets, respectively. The TIMN model has been extended from an earlier version treating binary IMN (BIMN) by adding NH₃ into the nucleation system and using a discrete-sectional bin structure to represent the sizes of clusters–particles starting from a single molecule up to background particles larger than a few micrometers.

2.2 Model representation of kinetic ternary nucleation processes

Figure 1 schematically illustrates the evolution of charged and neutral clusters–droplets explicitly simulated in the kinetic H₂SO₄–H₂O–NH₃ TIMN model. Here, H₂SO₄ (S) is the key atmospheric nucleation precursor driving the TIMN process while ions, H₂O (W), and NH₃ (A) stabilize the H₂SO₄ clusters and enhance H₂SO₄ nucleation rates in this way. Ions also enhance cluster formation rates due to the interaction with polar-nucleating species, leading to enhanced collision cross sections (Nadykto and Yu, 2003). The airborne ions are generated by galactic cosmic rays (GCRs) or produced by radioactive emanations, lightning, corona discharge, combustion, and other ionization sources. The initial negative ions, which are normally assumed to be ${\mathrm{NO}}_{\mathrm{3}}^{-}$, are converted into ${\mathrm{HSO}}_{\mathrm{4}}^{-}$ core ions (i.e., S⁻) and then to larger H₂SO₄ clusters in the presence of gaseous H₂SO₄. The initial positive ions H⁺W_w are converted into ${{\mathrm{H}}^{+}{\mathrm{A}}_{\mathrm{1}-\mathrm{2}}\mathrm{W}}_{\mathrm{w}}$ in the presence of NH₃, H⁺S_s W_w in the presence of H₂SO₄, or H⁺A_a S_sW_w in the case that both NH₃ and H₂SO₄ are present in the nucleating vapors. Some of the binary H₂SO₄–H₂O clusters, both neutral and charged, transform into ternary ones by taking up NH₃ vapors. The molar fraction of ternary clusters in nucleating vapors depends on [NH₃], the binding strength of NH₃ to binary and ternary pre-nucleation clusters, cluster composition, and ambient conditions such as T and RH.

Similar to the kinetic binary IMN (BIMN) model (Yu, 2006b), the kinetic TIMN model employs a discrete-sectional bin structure to represent clusters–particles. The bin index i represents the amount of core component (i.e., H₂SO₄) and i_d is the number of discrete bins. For small clusters (i≤i_d = 30 in this study), i is the number of H₂SO₄ molecules in the cluster (i.e., i=s) and the core volume of the ith bin $v_i=i\times v_{\mathrm{1}}$, where v₁ is the volume of one H₂SO₄ molecule. When i > i_d, v_i = VRAT_i × v_i−1, where VRAT_i is the volume ratio of the ith bin to the (i−1)th bin. The discrete-sectional bin structure enables the model to cover a wide range of sizes of nucleating clusters–particles with the highest possible size resolution for small clusters (Yu, 2006b). For clusters with a given bin i, the associated amounts of water and NH₃ and thus the effective radius of each ternary cluster are calculated based on the equilibrium of clusters–particles with the water vapor and/or ammonia, as described in later sections.

The evolution of positive, negative, and neutral clusters due to the simultaneous condensation, evaporation, recombination, coagulation, and other loss processes is described by the following differential equations obtained by the modification of those describing the evolution of the binary H₂SO₄–H₂O system (Yu, 2006b):

$$\begin{array}{ll}\text{(1)}& {\displaystyle }& {\displaystyle \frac{\partial N_{\mathrm{0}}^{+}}{\partial t}}=Q+{\mathit{\gamma }}_{\mathrm{1}}^{+}{N}_{\mathrm{1}}^{+}-N_{\mathrm{0}}^{+}{\displaystyle }& {\displaystyle }\cdot \left(\sum _{j=\mathrm{1}}^{i_{max}}{\mathit{\beta }}_{i,j}^{+}{N}_j^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i_{max}}{\mathit{\eta }}_{i,j}^{+}{N}_j^{+}+\sum _{j=\mathrm{0}}^{i_{max}}{\mathit{\alpha }}_{\mathrm{0},j}^{+,-}{N}_j^{-}\right)-N_{\mathrm{0}}^{+}{L}_{\mathrm{0}}^{+},\\ \text{(2)}& {\displaystyle }& {\displaystyle \frac{\partial N_{\mathrm{0}}^{-}}{\partial t}}=Q+{\mathit{\gamma }}_{\mathrm{1}}^{-}{N}_{\mathrm{1}}^{-}-N_{\mathrm{0}}^{-}{\displaystyle }& {\displaystyle }\cdot \left(\sum _{j=\mathrm{1}}^{i_{max}}{\mathit{\beta }}_{i,j}^{-}{N}_j^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i_{max}}{\mathit{\eta }}_{i,j}^{-}{N}_j^{-}+\sum _{j=\mathrm{0}}^{i_{max}}{\mathit{\alpha }}_{\mathrm{0},j}^{-,+}{N}_j^{+}\right)-N_{\mathrm{0}}^{-}{L}_{\mathrm{0}}^{-},\\ \text{(3)}& {\displaystyle }& {\displaystyle \frac{\partial N_{\mathrm{1}}^{\mathrm{0}}}{\partial t}}=P_{{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}}+\sum _{j=\mathrm{2}}^{i_{max}}{\mathit{\delta }}_{j,\mathrm{2}}{\mathit{\gamma }}_j^{\mathrm{0}}{N}_j^{\mathrm{0}}+\sum _{j=\mathrm{1}}^{i_{max}}\left({\mathit{\gamma }}_j^{+}{N}_j^{+}+{\mathit{\gamma }}_j^{-}{N}_j^{-}\right){\displaystyle }& {\displaystyle }-N_{\mathrm{1}}^{\mathrm{0}}\left(\sum _{j=\mathrm{1}}^{i_{max}}(\mathrm{1}-f_{\mathrm{1},j,\mathrm{1}}){\mathit{\beta }}_{\mathrm{1},j}^{\mathrm{0}}{N}_j^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i_{max}}\left({\mathit{\beta }}_{j,\mathrm{1}}^{+}{N}_j^{+}+{\mathit{\beta }}_{j,\mathrm{1}}^{-}{N}_j^{-}\right)\right)\\ {\displaystyle }& {\displaystyle }-N_{\mathrm{1}}^{\mathrm{0}}{L}_{\mathrm{1}}^{\mathrm{0}},\end{array}$$

$$\begin{array}{ll}\text{(4)}& {\displaystyle }& {\displaystyle \frac{\partial N_i^{+}\left(i\ge \mathrm{1}\right)}{\partial t}}=g_{i+\mathrm{1},i}{\mathit{\gamma }}_{i+\mathrm{1}}^{+}{N}_{i+\mathrm{1}}^{+}-g_{i,i-\mathrm{1}}{\mathit{\gamma }}_i^{+}{N}_i^{+}{\displaystyle }& {\displaystyle }+\sum _{j=\mathrm{0}}^{i-\mathrm{1}}\sum _{k=\mathrm{1}}^i{\displaystyle \frac{v_j}{v_i}}{f}_{j,k,i}{\mathit{\beta }}_{j,k}^{+}{N}_j^{+}{N}_k^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i-\mathrm{1}}\sum _{k=\mathrm{0}}^i{\displaystyle \frac{v_j}{v_i}}{f}_{j,k,i}{\mathit{\eta }}_{j,k}^{+}{N}_j^{+}{N}_k^{+}\\ {\displaystyle }& {\displaystyle }+\sum _{j=\mathrm{0}}^i\sum _{k=\mathrm{1}}^i{\displaystyle \frac{v_k}{v_i}}{f}_{j,k,i}{\mathit{\beta }}_{j,k}^{+}{N}_j^{+}{N}_k^{\mathrm{0}}-N_i^{+}\\ {\displaystyle }& {\displaystyle }\cdot \left(\sum _{j=\mathrm{1}}^{i_{max}}\left(\mathrm{1}-f_{i,j,i}\right){\mathit{\beta }}_{i,j}^{+}{N}_j^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i_{max}}\left(\mathrm{1}-f_{i,j,i}\right){\mathit{\eta }}_{i,j}^{+}{N}_j^{+}\right.\\ {\displaystyle }& {\displaystyle }\left.+\sum _{j=\mathrm{0}}^{i_{max}}{\mathit{\alpha }}_{i,j}^{+,-}{N}_j^{-}\right)-N_i^{+}{L}_i^{+},\end{array}$$

$$\begin{array}{ll}\text{(5)}& {\displaystyle }& {\displaystyle \frac{\partial N_i^{-}\left(i\ge \mathrm{1}\right)}{\partial t}}=g_{i+\mathrm{1},i}{\mathit{\gamma }}_{i+\mathrm{1}}^{-}{N}_{i+\mathrm{1}}^{-}-g_{i,i-\mathrm{1}}{\mathit{\gamma }}_i^{-}{N}_i^{-}{\displaystyle }& {\displaystyle }+\sum _{j=\mathrm{0}}^{i-\mathrm{1}}\sum _{k=\mathrm{1}}^i{\displaystyle \frac{v_j}{v_i}}{f}_{j,k,i}{\mathit{\beta }}_{j,k}^{-}{N}_j^{-}{N}_k^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i-\mathrm{1}}\sum _{k=\mathrm{0}}^i{\displaystyle \frac{v_j}{v_i}}{f}_{j,k,i}{\mathit{\eta }}_{j,k}^{-}{N}_j^{-}{N}_k^{-}\\ {\displaystyle }& {\displaystyle }+\sum _{j=\mathrm{0}}^i\sum _{k=\mathrm{1}}^i{\displaystyle \frac{v_k}{v_i}}{f}_{j,k,i}{\mathit{\beta }}_{j,k}^{-}{N}_j^{-}{N}_k^{\mathrm{0}}-N_i^{-}\left(\sum _{j=\mathrm{1}}^{i_{max}}\left(\mathrm{1}-f_{i,j,i}\right){\mathit{\beta }}_{i,j}^{-}{N}_j^{\mathrm{0}}\right.\\ {\displaystyle }& {\displaystyle }\left.+\sum _{j=\mathrm{0}}^{i_{max}}\left(\mathrm{1}-f_{i,j,i}\right){\mathit{\eta }}_{i,j}^{-}{N}_j^{-}+\sum _{j=\mathrm{0}}^{i_{max}}{\mathit{\alpha }}_{i,j}^{-,+}{N}_j^{+}\right)-N_i^{-}{L}_i^{-},\end{array}$$

$$\begin{array}{ll}\text{(6)}& {\displaystyle }& {\displaystyle \frac{\partial N_i^{\mathrm{0}}\left(i\ge \mathrm{2}\right)}{\partial t}}=g_{i+\mathrm{1},i}{\mathit{\gamma }}_{i+\mathrm{1}}^{\mathrm{0}}{N}_{i+\mathrm{1}}^{\mathrm{0}}-g_{i,i-\mathrm{1}}{\mathit{\gamma }}_i^{\mathrm{0}}{N}_i^{\mathrm{0}}{\displaystyle }& {\displaystyle }+\sum _{j=\mathrm{1}}^i\sum _{k=\mathrm{1}}^{i-\mathrm{1}}{\displaystyle \frac{v_k}{v_i}}{f}_{j,k,i}{\mathit{\beta }}_{j,k}^{\mathrm{0}}{N}_j^{\mathrm{0}}{N}_k^{\mathrm{0}}+\sum _{j=\mathrm{0}}^i\sum _{k=\mathrm{0}}^i{f}_{j,k,i}{\mathit{\alpha }}_{j,k}^{+,-}\\ {\displaystyle }& {\displaystyle }\cdot \left({\displaystyle \frac{v_k}{v_i}}{N}_j^{+}{N}_k^{-}+{\displaystyle \frac{v_j}{v_i}}{N}_j^{+}{N}_k^{-}\right)-N_i^{\mathrm{0}}\\ {\displaystyle }& {\displaystyle }\cdot \left(\sum _{j=\mathrm{1}}^{i_{max}}\left(\mathrm{1}-f_{i,j,i}\right){\mathit{\beta }}_{i,j}^{\mathrm{0}}{N}_j^{\mathrm{0}}+\sum _{j=\mathrm{0}}^{i_{max}}\left({\mathit{\beta }}_{j,i}^{+}{N}_j^{+}+{\mathit{\beta }}_{j,i}^{-}{N}_j^{-}\right)\right)\\ {\displaystyle }& {\displaystyle }-N_i^{\mathrm{0}}{L}_i^{\mathrm{0}}.\end{array}$$

In Eqs. (1)–(6), the superscripts “+”, “−”, and “0” refer to positive, negative, and neutral clusters, respectively, while subscripts i, j, and k represent the bin indexes. $N_{\mathrm{0}}^{+,-}$ and Q are the concentration of initial ions not containing H₂SO₄ (i.e., H⁺A_aW_w and ${\mathrm{NO}}_{\mathrm{3}}^{-}$) and the ionization rate, respectively. N_i is the total number concentration (cm⁻³) of all cluster/particles (binary + ternary) in bin i. For small clusters (i≤i_d), N_i is the number concentration (cm⁻³) of all clusters containing i H₂SO₄ molecules. For example, $N_{\mathrm{1}}^{\mathrm{0}}$ is the total concentration of binary and ternary neutral clusters containing one H₂SO₄ molecule. Index i in Eq. (4) refers to the sum of H₂SO₄ and ${\mathrm{HSO}}_{\mathrm{4}}^{-}$. The second term of Eq. (2) describes the reaction of ${\mathrm{HSO}}_{\mathrm{4}}^{-}$ + HNO₃ → ${\mathrm{NO}}_{\mathrm{3}}^{-}$ + H₂SO₄. Although the rate of this reaction is generally negligible, we keep the term there for completeness. $P_{{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}}$ is the gas-phase production rate of neutral H₂SO₄ molecules. $L_i^{+,-,\mathrm{0}}$ is the loss rate due to scavenging by preexisting particles and wall and dilution losses in the laboratory chamber studies (Kirkby et al., 2011; Olenius et al., 2013; Kürten et al., 2016). $f_{j,k,i}$ is the volume fraction of intermediate particles (volume = v_j + v_k) partitioned into bin i with respect to the core component – H₂SO₄, as defined in Jacobson et al. (1994). $g_{i+\mathrm{1},i}=v_{\mathrm{1}}/(v_{i+\mathrm{1}}-v_i)$ is the volume fraction of intermediate particles of volume ($v_{i+\mathrm{1}}-v_{\mathrm{1}}$) partitioned into bin i. ${\mathit{\delta }}_{j,\mathrm{2}}=\mathrm{2}$ at j=2 and ${\mathit{\delta }}_{j,\mathrm{2}}=\mathrm{1}$ at j≠2. ${\mathit{\gamma }}_i^{+}$, ${\mathit{\gamma }}_i^{-}$ , and ${\mathit{\gamma }}_i^{\mathrm{0}}$ are the mean (or effective) cluster evaporation coefficients for positive, negative and neutral clusters in bin i, respectively. ${\mathit{\beta }}_{i,j}^{+}$, ${\mathit{\beta }}_{i,j}^{-}$, and ${\mathit{\beta }}_{i,j}^{\mathrm{0}}$ are the coagulation kernels for the neutral clusters–particles in bin j interacting with positive, negative, and neutral clusters–particles in bin i, respectively, which reduce to the condensation coefficients for H₂SO₄ monomers at j=1. ${\mathit{\eta }}_{j,k}^{+}$ and ${\mathit{\eta }}_{j,k}^{-}$ are coagulation kernels for clusters–particles of like sign from bin j and clusters–particles from bin k. It should be noted that the electrostatic repulsion is too strong for small clusters to gain more than one charge. However, small charged clusters can be scavenged by large preexisting particles of the same polarity. Large preexisting particles serve as the sink for small clusters in the model and the effect of multiple charge is small and thus is not tracked. ${\mathit{\alpha }}_{i,j}^{+,-}$ is the recombination coefficient for positive clusters–particles in bin i interacting with negative clusters–particles in bin j, while ${\mathit{\alpha }}_{i,j}^{-,+}$ is the recombination coefficient of negative clusters–particles from bin i interacting with positively charged clusters–particles from bin j.

The methods for calculating β, γ, η, and α for binary H₂SO₄–H₂O clusters have been described in our previous publications (Yu and Turco, 2001; Nadykto and Yu, 2003; Yu, 2006b). Dipole–charge interaction (Nadykto and Yu, 2003), image capture, and three-body trapping effects (Hoppel and Frick, 1986) are considered in the calculation of these coefficients. Since β, η, and α depend on the cluster mass (or size) rather than on the cluster composition, schemes for calculating these properties in binary and ternary clusters are identical. In contrast, γ is quite sensitive to cluster composition. The evaporation rate coefficient of H₂SO₄ molecules from clusters containing i H₂SO₄ molecules (γ_i) is largely controlled by the stepwise Gibbs free energy change $\mathrm{\Delta }{G}_{i-\mathrm{1},i}$ of formation of an i-mer from an (i−1)-mer (Yu, 2007)

$$\begin{array}{}\text{(7)}& {\displaystyle }& {\displaystyle }{\mathit{\gamma }}_i={\mathit{\beta }}_{i-\mathrm{1}}{N}^o\exp \left({\displaystyle \frac{\mathrm{\Delta }{G}_{i-\mathrm{1},i}}{RT}}\right),\text{(8)}& {\displaystyle }& {\displaystyle }\mathrm{\Delta }{G}_{i-\mathrm{1},i}=\mathrm{\Delta }{H}_{i-\mathrm{1},i}^o-T\mathrm{\Delta }{S}_{i-\mathrm{1},i}^o,\end{array}$$

where R is the molar gas constant, N^o is the arbitrary number concentration of a hypothetical gas consisting solely of the species for which the calculation is performed (generally under the reference vapor pressure P of 1 atm). ΔH^o and ΔS^o are enthalpy and entropy changes under the standard conditions (T=298 K, P=1 atm), respectively. The temperature dependence of ΔH^o and ΔS^o, which is generally small and typically negligible over the temperature range of interest (Nadykto et al., 2009), was not considered.

Figure 2 Stepwise Gibbs free energy change under standard conditions for the addition of a water ($\mathrm{\Delta }{G}_{+\mathrm{W}}^o$), ammonia ($\mathrm{\Delta }{G}_{+\mathrm{A}}^o$), or sulfuric acid ($\mathrm{\Delta }{G}_{+\mathrm{S}}^o$) molecule to form the given positively charged clusters as a function of the number of water molecules in the clusters (w). Lines are QC-based values, and symbols are experimental results or semi-experimental estimates (see notes under Table A2 for the references).

2.3 Thermochemical data of neutral and charged binary and ternary clusters

ΔH, ΔS, and ΔG values needed to calculate cluster evaporation rates (Eq. 7) for the TIMN model can be derived from laboratory measurements and computational quantum chemistry (QC) calculation. Thermochemical properties of neutral and charged binary and ternary clusters obtained using the computational chemical methods and comparisons of computed energies with available experimental data and semi-experimental estimates are given in Tables A1–A4 and discussed in the Appendix. As an example, Fig. 2 shows ΔG associated with the addition of water ($\mathrm{\Delta }{G}_{+\mathrm{W}}^o$), ammonia ($\mathrm{\Delta }{G}_{+\mathrm{A}}^o$), and sulfuric acid ($\mathrm{\Delta }{G}_{+\mathrm{S}}^o$) to binary and ternary clusters as a function of the cluster hydration number w. H₂O has high proton affinity and thus H₂O is strongly bonded to all positive ions with low w. $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ expectedly becomes less negative and binding of H₂O to binary and ternary clusters weakens due to the screening effect as the hydration number w grows (Fig. 2a). The presence of NH₃ in the clusters weakens binding of H₂O to positive ions. For example, $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ for H⁺A₁W_wS₁ is ∼3–4 kcal mol⁻¹ less negative than that for H⁺W_wS₁ at w=3–6. The addition of one more NH₃ to the clusters to form H⁺A₂W_w and H⁺A₂W_wS₁ further weakens H₂O binding by ∼1.5–6 kcal mol⁻¹ at w=1–3, while exhibiting a much smaller impact on hydration free energies at w > 3. Both the absolute values and trends in $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ derived from calculations are in agreement with the laboratory measurements within the uncertainty range of ∼1–2 kcal mol⁻¹ for both QC calculations and measurements. This confirms the efficiency and precision of QC methods in calculating thermodynamic data needed for the development of nucleation models.

The proton affinity of NH₃ is 204.1 kcal mol⁻¹, which is 37.5 kcal mol⁻¹ higher than that of H₂O (166.6 kcal mol⁻¹) (Jolly, 1991). The hydrated hydronium ions (H⁺W_w) are easily converted to H⁺A₁W_w in the presence of NH₃. The binding of NH₃ and H₂O molecules to H⁺W_w exhibits a similar pattern. In particular, binding of NH₃ to H⁺W_w decreases as w is growing, with $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ for H⁺A₁W_w ranging from −52.08 kcal mol⁻¹ at w=1 to −8.32 kcal mol⁻¹ at w=9. The binding of NH₃ to H⁺W_wS₁ ions is also quite strong, with $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ for H⁺A₁W_wS₁ ranging from −33.14 kcal mol⁻¹ at w=1 and to −10.57 kcal mol⁻¹ at w=6. The addition of the NH₃ molecule to H⁺A₁W_w (to form H⁺A₂W_w) is much less favorable thermodynamically than that to H⁺W_w, with the corresponding $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ being −22 and −6 kcal mol⁻¹ at w=2 and w=6, respectively. The $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ values for H⁺A₂W_w are 3–5 kcal mol⁻¹ more negative than the experimental values at w=0–1; however, they are pretty close to experimental data at w=2–3 (Fig. 2b and Table A2). While it is possible that the QC method overestimates the charge effect on the formation free energies of smallest clusters, the possible overestimation at w=0–1 will not affect nucleation calculations because most H⁺A₂W_w clusters in the atmosphere contain more than two water molecules (i.e., w > 2) due to the strong hydration (see Table A2 and Fig. 2a).

A comparison of QC and semi-experimental estimates of $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ values associated with the attachment of H₂SO₄ to positive ions shown in Fig. 2c indicates that computed $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ values agree well with observations for H⁺W_wS₁ and H⁺A₁W_wS₁ but differ by ∼2–4 kcal mol⁻¹ from semi-experimental values for H⁺A₂W_wS₁. As seen from Fig. 2a and c, the attachment of NH₃ to H⁺W_wS₁ weakens the binding of both H₂O and H₂SO₄ to the clusters. This suggests that the attachment of NH₃ leads to the evaporation of H₂SO₄ and H₂O molecules from the clusters. In other words, H₂SO₄ is less stable in H⁺A₁W_wS₁ than in H⁺W_wS₁ (Fig. 2c). While this may be taken for the indication that NH₃ inhibits nucleation on positive ions at the first look, further calculations show that binding of NH₃ to H⁺A₁W_wS₁ is quite strong (Fig. 2b) and that H₂SO₄ in the H⁺A₂W_wS₁ cluster is much more stable than that in H⁺A₁W_wS₁, with $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ being more negative by ∼7 kcal mol⁻¹ at w > 2. The H⁺A₂W_wS₁ cluster can also be formed via the attachment of H₂SO₄ to H⁺A₂W_w. In the presence of sufficient concentrations of NH₃, a large fraction of positively charged H₂SO₄ monomers exist in the form of H⁺A₂W_wS₁ and hence NH₃ enhances nucleation of positive ions. Since positively charged H₂SO₄ dimers are expected to contain a large number of water molecules, we have not yet computed and derived quantum chemical data for these clusters. The CLOUD measurements do indicate that once H⁺A₂W_wS₁ clusters are formed, they can continue to grow to larger H⁺A_aW_wS_s clusters along $a=s+\mathrm{1}$ pathway (Schobesberger et al., 2015).

Figure 2 clearly shows that the calculated values in most cases agree with measurements within the uncertainty range that justifies the application of QC values in the case that no reliable experimental data are available.

2.4 Nucleation barriers for neutral and charged clusters and size-dependent evaporation rates

Nucleation barriers and cluster evaporation rates are critically important for calculations of nucleation rates. This section describes the methods employed to calculate the evaporation rates of nucleating clusters of variable sizes and compositions (i.e., γ in Eqs. 1–6) in the TIMN model.

2.4.1 Equilibrium distributions of small binary and ternary clusters

In the atmosphere, [H₂O] is much higher than [H₂SO₄] and thus H₂SO₄ clusters–particles are always in equilibrium with water vapor (Yu, 2007). In the lower troposphere, where most of the nucleation events were observed, [H₂SO₄] is typically at a sub-parts per trillion to parts per trillion level, while [NH₃] is in the range of sub-parts per billion to parts per billion levels (Butler et al., 2016; Warner et al., 2016) (note that, in what follows, all references to vapor mixing ratios – parts per billion and parts per trillion – are by volume). This means that small ternary clusters can be considered to be in equilibrium with H₂O and NH₃ vapors. Like the previous BIMN model derived assuming equilibrium of binary clusters with water vapor, the present TIMN model treats small clusters containing a given number of H₂SO₄ molecules as being in equilibrium with both H₂O and NH₃. Their relative concentrations are calculated using the thermodynamic data shown in Tables A1–A4. It should be noted that the system may deviate from equilibrium and the model scheme is probably not suitable when [NH₃] is less than or close to [H₂SO₄]. Under such cases, the equilibrium assumption may overestimate nucleation rates.

Figure 3 Relative abundance (or molar fraction) of small clusters containing a given number of H₂SO₄ molecules for positive (a), negative (b), and neutral (c, d) cluster types at a temperature of 292 K and three different combinations of RHs (38 % and 90 %) and [NH₃] (0.3 and 3 ppb). Some clusters with different numbers of water molecules were grouped together to make the plot more clear and neat. For the clusters shown in (d), there are no hydrate data and thus hydration for these clusters was not calculated.

Figure 3 shows the relative abundance (or molar fractions) of small positive, negative, and neutral clusters ($f_{s,a,w}^{+,-,\mathrm{0}}$) containing a given number of H₂SO₄ molecules at the ambient temperature of 292 K and three different combinations of RH and [NH₃] values. As a result of relative instability of H₂SO₄ in H⁺A₁W_wS₁ compared to H⁺W_wS₁ or H⁺A₂W_wS₁ (Fig. 2c), most of positive ions with one H₂SO₄ molecule exist in the form of either H⁺W_wS₁ or H⁺A₂W_wS₁ (i.e, containing either zero or two NH₃ molecules; Fig. 3a). When [NH₃] = 0.3 ppb (with T=292 K), most of the positive ions containing one H₂SO₄ molecule do not contain NH₃ and their composition is dominated by H⁺W_wS₁ ($\stackrel{\mathrm{‾}}{w}=\sim \mathrm{7}$). At the given T and [NH₃] = 0.3 ppb, around 17 % of positive ions with one H₂SO₄ molecule contain two NH₃ molecules at RH = 38 %. The fraction of positive ions containing one H₂SO₄ and two NH₃ molecules decreases to 0.9 %, when RH=90 %. At T=292 K and RH=38 %, the increase in [NH₃] by a factor of 10 to 3 ppb leads to the domination of H⁺A₂W_wS₁ (∼95 %) in the composition of positively charged H₂SO₄ monomers. As expected, the composition of positive ions and their contribution to nucleation depends on T, RH, and [NH₃]. The incorporation of the quantum chemical and experimental clustering thermodynamics in the framework of the kinetic nucleation model enables us to study all these dependencies.

As a result of very weak binding of H₂O and NH₃ to small negative ions (Table A4), nearly all negatively charged clusters with s = 0–1 do not contain water and ammonia (not shown). In the case that s grows to 2, all S⁻S₂A_aW_w clusters still do not contain NH₃ (i.e., a=0), while only 20 %–40 % of them contain one water molecule (w=1) (Fig. 3b). As s further increases to 3, NH₃ begins to enter some of the negatively charged ions. The fraction of S⁻S₃A_aW_w clusters containing one NH₃ molecule is 9 % at RH = 38 % and [NH₃] = 0.3 ppb, 3 % at RH = 90 % and [NH₃] = 0.3 ppb, and 50 % at RH = 38 % and [NH₃] = 3 ppb. Most S⁻S₃W_w clusters are hydrated while the fraction of S⁻S₃A_aW_w clusters containing two NH₃ molecules at these ambient conditions is negligible. The fraction of negative cluster ions containing two NH₃ molecules becomes significant at s=4 (Fig. 3b) and increases from 28 % at [NH₃] = 0.3 ppb to 80 % at [NH₃] = 3 ppb at RH = 38 %. At [NH₃] = 0.3 ppb, the increase in RH from 38 % to 90 % reduces the fraction of NH₃-containing S⁻S₃A_aW_w clusters (i.e, a > =1) from to 95 % to 70 %, demonstrating a significant impact of RH on cluster compositions and emphasizing the importance of accounting for the RH in calculations of ternary nucleation rates.

The equilibrium distributions of neutral clusters are presented in Fig. 3c (H₂SO₄ monomers and dimers) and Fig. 3d (H₂SO₄ trimers and tetramers). Hydration is accounted for in the case of monomers and dimers and not included, due to lack of thermodynamic data, in calculations for trimers and tetramers. Based on the thermodynamic data shown in Table A3, the dominant fraction of neutral monomers is hydrated (79 % at RH = 38 % and 94 % at RH = 90 %) while the fraction of monomers containing NH₃ is negligible (0.02 % at [NH₃] = 0.3 ppb and 0.2 % at [NH₃] = 3 ppb, RH = 38 %). As a result of the growing binding strength of NH₃ with the cluster size (Table A3), the fraction of neutral sulfuric acid dimers containing one NH₃ molecule reaches 18 % at [NH₃] = 0.3 ppb and 69 % at [NH₃] = 3 ppb when T = 292 K and RH = 38 %. In the case of H₂SO₄ trimers and tetramers, data shown in Fig. 3d are limited to the relative abundance of unhydrated clusters only. Under the given conditions, most trimers contain two NH₃ molecules while most tetramers contain three NH₃ molecules. At [NH₃] = 3 ppb, ∼2 % of trimers contain three NH₃ molecules (i.e., $s=a=\mathrm{3}$) and 55 % of tetramers contain four NH₃ molecules (i.e., $s=a=\mathrm{4}$). As a result of a significant drop of $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ in the case that the a∕s ratio exceeds one (Table A3), the fraction of neutral clusters with $a=s+\mathrm{1}$ is negligible. The cluster distributions clearly indicate that small sulfuric acid clusters are still not fully neutralized by NH₃ even if [NH₃] is at a parts per billion level and that the degree of neutralization (i.e., a:s ratio) increases with the cluster size.

2.4.2 Mean stepwise and accumulative Gibbs free energy change and impact of ammonia

In the TIMN model, the equilibrium distributions are used to calculate number-concentration-weighted stepwise Gibbs free energy change for adding one H₂SO₄ molecule to form a neutral, positively charged, and negatively charged cluster containing s H₂SO₄ molecules (${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$):

$$\begin{array}{}\text{(9)}& {\displaystyle }{\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}^{+,-,\mathrm{0}}={\sum }_{a,w}{f}_{s,a,w}^{+,-,\mathrm{0}}\mathrm{\Delta }{G}_{s-\mathrm{1},s,a,w}^{+,-,\mathrm{0}},\end{array}$$

where $f_{s,a,w}^{+,-,\mathrm{0}}$ is the equilibrium fraction of a particular cluster within a cluster type as shown in Fig. 3.

In the atmosphere, where substantial nucleation is observed, the sizes of critical clusters are generally small (s < ∼5–10) (e.g., Sipilä et al., 2010) and nucleation rates are largely controlled by the stability (or γ) of small clusters with s < ∼5–10. QC calculations and experimental data on clustering thermodynamics available for clusters of small sizes (Tables A2–A4) are critically important as the formation of these small clusters is generally the limiting step for nucleation. Nevertheless, thermodynamics data for larger clusters are also needed to develop a robust nucleation model that can calculate nucleation rates under various conditions. Both measurements and QC calculations (Tables A2–A4) show significant effects of charge and charge signs (i.e., positive or negative) on the stability and composition of small clusters. These charge effects decrease quickly as the clusters grow due to the short-ranged nature of dipole–charge interaction and the quick decrease in electrical field strength around charged clusters as cluster sizes increase (Yu, 2005). Based on experimental data (Kebarle et al., 1967; Davidson et al., 1977; Wlodek et al., 1980; Holland and Castleman, 1982; Froyd and Lovejoy, 2003b), the stepwise ΔG values for clusters decrease exponentially as the cluster sizes increase and approach to the bulk values when clusters contain more than ∼8–10 molecules (Yu, 2005). Cluster compositions measured with an atmospheric pressure interface time-of-flight (APi-TOF) mass spectrometer during CLOUD experiments also show that the difference in the composition of positively and negatively charged clusters quickly decreases as the number of H₂SO₄ molecules increases from 1 to ∼10 and exhibits little further changes (Schobesberger et al., 2015).

In the present TIMN model, we assume that both neutral and charged clusters have the same composition when s≥10 and the following extrapolation scheme is used to calculate $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ for clusters up to s=10:

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }\mathrm{\Delta }{G}_{s-\mathrm{1},s}=\mathrm{\Delta }{G}_{s_{\mathrm{1}}-\mathrm{1},s_{\mathrm{1}}}\\ \text{(10)}& {\displaystyle }& {\displaystyle }+{\displaystyle \frac{\left(\mathrm{\Delta }{G}_{s_{\mathrm{2}}-\mathrm{1},s_{\mathrm{2}}}-\mathrm{\Delta }{G}_{s_{\mathrm{1}}-\mathrm{1},s_{\mathrm{1}}}\right)\left(e^{-sc}-e^{-s_{\mathrm{1}}c}\right)}{\left(e^{-s_{\mathrm{2}}c}-e^{-s_{\mathrm{1}}c}\right)}},\end{array}$$

where $\mathrm{\Delta }{G}_{s_{\mathrm{1}}-\mathrm{1},s_{\mathrm{1}}}$ is the stepwise mean Gibbs free energy change for H₂SO₄ addition for a specific type (neutral, positive, or negative) of clusters at s=s₁ that can be derived from QC calculation and/or experimental measurements, and $\mathrm{\Delta }{G}_{s_{\mathrm{2}}-\mathrm{1},s_{\mathrm{2}}}$ is the corresponding value for clusters at s=s₂ (= 10 in the present study) that is calculated in the capillarity approximation accounting for the Kelvin effect. c in Eq. (10) is the exponential coefficient that determines how fast $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ approaches to bulk values as s increases. In the present study, c is estimated by fitting $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ at s=2 and s=3 based on Eq. (10) to the corresponding $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ from experimental (Hanson and Lovejoy, 2006; Kazil et al., 2007) or QC data (Table A3). Apparently the interpolation approximation Eq. (10) is subject to uncertainty. Nevertheless, it is a reasonable approach to connect thermochemical properties of QC data for small binary and ternary clusters that cannot be adequately described by the capillarity approximation with those for large clusters that can be adequately described by the very same capillarity approximation, and it is the best approach we can come up with at this point in order to develop a model that can be applied to all conditions. Further QC and experimental studies of the thermodynamics of relatively larger clusters can help to reduce the uncertainty.

For clusters with s≥s₂, the capillarity approximation is used to calculate $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ as

$$\begin{array}{}\text{(11)}& {\displaystyle }\mathrm{\Delta }{G}_{s-\mathrm{1},s}=-RT\ln (P/P_{\mathrm{s}})+{\displaystyle \frac{\mathrm{2}\mathit{\sigma }{v}_{\mathrm{1}}{N}_{\mathrm{A}}}{r_{\mathrm{s}}}},\end{array}$$

where P is the H₂SO₄ vapor pressure and P_s is the H₂SO₄ saturation vapor pressure over a flat surface with the same composition as the cluster. σ is the surface tension and v₁ is the volume of one H₂SO₄ molecule. r_s is the radius of the cluster and N_A is the Avogadro's number.

The scheme to calculate bulk $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ (s≥10) for H₂SO₄–H₂O binary clusters has been described in Yu (2007). For ternary nucleation, both experiments (Schobesberger et al., 2015) and QC calculations (Table A4) indicate that the growth of relatively large clusters follows the s=a line (i.e, in the composition of ammonia bisulfate). In the present TIMN model, the bulk $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ values for ternary clusters are calculated based on parameterized H₂SO₄ saturation vapor pressure over ammonia bisulfate as a function of temperature, derived by Marti et al. (1997) from vapor pressures measured at a temperature between 27 and 60 ^∘C and surface tension measured at 298 K from Hyvärinen et al. (2005). The uncertainty in saturation vapor pressures and surface tension used in the calculation of the bulk $\mathrm{\Delta }{G}_{s-\mathrm{1},s}$ values is another source of uncertainty in the TIMN model, although it is likely to be small compared to other uncertainties as the nucleation is generally limited by the formation of small clusters.

Figure 4 (a) Average stepwise Gibbs free energy change for the addition of one H₂SO₄ molecule to form a neutral (black), positively charged (red), or negatively charged (blue) binary H₂SO₄–H₂O (dashed lines or empty circles) or ternary H₂SO₄–H₂O–NH₃ (solid lines or filled circles) cluster containing s H₂SO₄ molecules (${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$). (b) Same as (a) but for the cumulative (total) Gibbs free energy change in each case. Filled and empty circles in (a) refer to ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ obtained using measurements and/or quantum-chemical calculations. ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ for larger clusters with s ≥ 10, which approach the properties of the equivalent bulk liquid (20), is calculated using the capillarity approximation. Interpolation is used to calculate ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ for clusters up to s=10 (Eq. 11). Calculations were carried out at T = 292 K, RH = 38 %, [H₂SO₄] = 3 × 10⁸ cm⁻³, and [NH₃] = 0.3 ppb. The inset diagrams represent equilibrium geometries for the most stable isomers of selected binary clusters ((H₃O⁺)(H₂SO₄)(H₂O)₆, (H₂SO₄)₂(H₂O)₄, and (${\mathrm{HSO}}_{\mathrm{4}}^{-}$)(H₂SO₄)₄(H₂O)₂) and ternary clusters ((${\mathrm{NH}}_{\mathrm{4}}^{+}$)(H₂SO₄)(NH₃)(H₂O)₄, (${\mathrm{HSO}}_{\mathrm{4}}^{-}$)(H₂SO₄)₄(H₂O)(NH₃), (H₂SO₄)₄(NH₃)₄).

Figure 4 presents stepwise (${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$) and cumulative (total) ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{\mathrm{s}}$ Gibbs free energy changes associated with the formation of neutral, positively charged, and negatively charged binary and ternary clusters containing s H₂SO₄ molecules under the conditions specified in the figure caption. The clusters are assumed to be in equilibrium with water (Yu, 2007) and ammonia (Fig. 3). As seen from Fig. 4, the presence of NH₃ reduces the mean ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ for larger clusters, which can be treated as the bulk binary H₂SO₄–H₂O solution (Schobesberger et al., 2015), by ∼3 kcal mol⁻¹, indicating a substantial reduction in the H₂SO₄ vapor pressure over ternary solutions (Marti et al., 1997). The comparison also shows that the influence of NH₃ on ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ of small clusters ($s\le \sim \mathrm{4}$) is much lower than that on larger ones and bulk solutions. For example, at [NH₃] = 0.3 ppb, the differences in ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ between binary and ternary positive ions with s=1 and neutral clusters with s=2 are only 0.45 and ∼1 kcal mol⁻¹, respectively. In the case of negative ions, zero and 0.27–0.45 kcal mol⁻¹ differences at s≤2 and s=3–4, respectively, were observed. The reduced effect of ammonia on smaller clusters is explained (Tables A2–A4) by ammonia's weaker bonding to smaller clusters than to larger ones, which in turn yields lower average NH₃-to-H₂SO₄ ratios (Fig. 3). It should be noted that QC data for positively charged clusters are very limited and the interpolation approximation is subject to large uncertainty. In order for the nucleation on positive ions to occur, the first step is for H₂SO₄ to attach to a positive ion that does not contain H₂SO₄. Unlike negative ions, the effect of charge on the bonding of H₂SO₄ with positive ions is much weaker and thus the stepwise Gibbs free energy change for the addition of one H₂SO₄ molecule to form a positively charged cluster is likely to be similar to that of neutral clusters, i.e., decreasing with cluster size. Therefore, the QC data for positively charged clusters containing one H₂SO₄ molecule provide a critical constraint. The success of the model in predicting the [NH₃] needed for nucleation on positive ions to occur (see Sect. 3) shows the usefulness of the first-step data and approximation.

As seen from Fig. 4, bonding of H₂SO₄ to small negatively charged clusters (s < 3) is much stronger than that to neutrals and positive ions. As a result, at s < 3 the formation of negatively charged clusters is barrierless (${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ < 0). These small clusters cannot be considered nucleated particles because ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ (Fig. 4a) first increases and then decreases with growing s, reaching the maximum barrier values at $s=\sim \mathrm{3}$–6. ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ can become positive for larger clusters due to the charge effect decreasing quickly as the clusters are growing. The effect of NH₃ on negative ions becomes important at $s\ge \sim \mathrm{4}$, when bonding between the clusters and NH₃ becomes strong enough to contaminate a large fraction of binary clusters with ammonia (Fig. 3). In contrast, the impact of NH₃ on neutral dimers and positively charged monomers of H₂SO₄, as well as on ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ for both positively charged and neutral clusters, monotonically decreases for all s values, including s≤5.

${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ for charged and neutral clusters converges into the bulk values at $s=\sim \mathrm{10}$, when impact of the chemical identity of the core ion on the cluster composition becomes diffuse (Schobesberger et al., 2015) and when the contribution of the electrostatic effect to ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ becomes less than ∼0.5 kcal mol⁻¹. The comparison of cumulative (total) ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{\mathrm{s}}$ (Fig. 4b) indicates the lowest nucleation barrier for the case of negative ions, followed by positive ions and neutrals. The barrierless formation of clusters with s ranging from 1 to 3 substantially reduces the nucleation barrier for negatively charged ions and facilitates their nucleation. The presence of 0.3 ppb of NH₃ lowers the nucleation barrier for negative, positive, and neutral clusters from ∼17, 24, and 38 kcal mol⁻¹ to 2, 7, and 16 kcal mol⁻¹, respectively. A relatively low nucleation barrier for charged ternary clusters is explained by the simultaneous effect of ionization and NH₃, which also reduces the size of the critical cluster (s^*).

It is important to note that the size of the critical cluster, commonly used to “measure” the activity of nucleation agents in the classical nucleation theory (Coffman and Hegg, 1995; Korhonen et al., 1999; Vehkamäki et al., 2002; Napari et al., 2002; Hamill et al., 1982) is no longer a valid indicator, when charged molecular clusters and small nanoparticles are considered. As seen from Fig. 4, positively charged ternary critical clusters ($s^{*}=\mathrm{3}$–4) are smaller than the corresponding negatively charged ones ($s^{*}=\mathrm{4}$–5); however, the nucleation barrier for ternary positive clusters under the condition specified in the figure caption is more than 3 times higher than that for ternary negatives ones.

2.4.3 Size- and composition-dependent H₂SO₄ evaporation rates

As we mentioned earlier, H₂SO₄ is the key atmospheric nucleation precursor driving the formation and growth of clusters in the ternary H₂SO₄–H₂O–NH₃ system while ions, H₂O, and NH₃ act to stabilize the H₂SO₄ clusters. The clustering thermodynamic data derived from QC calculations and measurements (Sect. 2.3) are used to constrain size- and composition- dependent Gibbs free energy changes and evaporation rates of H₂SO₄, which are critically important. Average or effective rates of H₂SO₄ molecule evaporation from positively charged, negatively charged, and neutral clusters containing s H₂SO₄ molecules (${\stackrel{\mathrm{‾}}{\mathit{\gamma }}}_{\mathrm{s}}^{+,-,\mathrm{0}}$) are calculated from ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ as

$$\begin{array}{}\text{(12)}& {\displaystyle }{\stackrel{\mathrm{‾}}{\mathit{\gamma }}}_{\mathrm{s}}^{+,-,\mathrm{0}}={\mathit{\beta }}_{s-\mathrm{1}}^{+,-,\mathrm{0}}{N}^o\exp \left({\displaystyle \frac{{\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}}{RT}}\right),\end{array}$$

where N^o is as defined in Eq. (6). The present model assumes only a single H₂SO₄ molecule evaporates, i.e., no water ligands, for instance, are attached to it. This is likely the dominant evaporation pathway as hydrated H₂SO₄ molecules are generally more stable.

Figure 5 The number-concentration-weighted mean evaporation rates ($\stackrel{\mathrm{‾}}{\mathit{\gamma }}$) of H₂SO₄ molecules from neutral clusters (black), positively charged clusters (red), and negatively charged clusters (blue) for binary (H₂SO₄–H₂O, dashed lines) and ternary (H₂SO₄–H₂O–NH₃, solid lines) nucleating systems containing s H₂SO₄ molecules (${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$). T = 292 K, RH = 38 %, and [NH₃] = 0.3 ppb for the ternary system.

Figure 5 gives the mean evaporation rate ($\stackrel{\mathrm{‾}}{\mathit{\gamma }}$) of a H₂SO₄ molecule from these clusters under the conditions corresponding to Fig. 4. The shapes of $\stackrel{\mathrm{‾}}{\mathit{\gamma }}$ curves are similar to those of ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ (Fig. 4a) as $\stackrel{\mathrm{‾}}{\mathit{\gamma }}$ values are largely controlled by ${\stackrel{\mathrm{‾}}{\mathrm{\Delta }G}}_{s-\mathrm{1},s}$ (Eq. 12). The presence of ammonia, as expected, significantly reduces the vapor pressure of H₂SO₄ over bulk aerosol (Marti et al., 1997) and hence the H₂SO₄ evaporation rate. The evaporation rates of both neutral and positive clusters decrease as s increases, and the positive clusters are uniformly more stable than corresponding neutral clusters. $\stackrel{\mathrm{‾}}{\mathit{\gamma }}$ for negative ions first increases and then decreases as s increases, peaking at around $s=\sim \mathrm{3}$–6. The presence of NH₃ reduces the evaporation rates of larger clusters by more than 2 orders of magnitude and the effect decreases for smaller clusters, as the binding of NH₃ to small neutral and charged clusters is weaker compared to that for larger clusters (Fig. 4). [NH₃] influences the average NH₃ : H₂SO₄ ratio (Fig. 3) and the evaporation rates of these small clusters. The nucleation rates, limited by formation of small clusters (s < ∼5), depend strongly on the stability or evaporation rate of these small clusters. While the binding of NH₃ to small neutral and charged clusters is weaker compared to that to larger clusters, small clusters containing NH₃ are much more stable than those without (Fig. 4) and thus ammonia is important for nucleation.

3 TIMN rates and comparisons with CLOUD measurements

The evolution of cluster/particle size distributions can be obtained by solving the dynamic Eqs. (1)–(6). Since the concentrations of clusters of all sizes are predicted, the nucleation rates in the kinetic model can be calculated for any cluster size larger than the critical size of neutral clusters (i > i^*) (Yu, 2006b),

$$\begin{array}{ll}{\displaystyle }{J}_i& {\displaystyle }=J_i^{+}+J_i^{-}+J_i^{\mathrm{0}}={\mathit{\beta }}_{i,\mathrm{1}}^{+}{N}_{\mathrm{1}}^{\mathrm{0}}{N}_i^{+}-{\mathit{\gamma }}_i^{+}{N}_{i+\mathrm{1}}^{+}+{\mathit{\beta }}_{i,\mathrm{1}}^{-}{N}_{\mathrm{1}}^{\mathrm{0}}{N}_i^{-}\\ \text{(13)}& {\displaystyle }& {\displaystyle }-{\mathit{\gamma }}_i^{-}{N}_{i+\mathrm{1}}^{-}+{\mathit{\beta }}_{i,\mathrm{1}}^{\mathrm{0}}{N}_{\mathrm{1}}^{\mathrm{0}}{N}_i^{\mathrm{0}}-{\mathit{\gamma }}_i^{\mathrm{0}}{N}_{i+\mathrm{1}}^{\mathrm{0}},\end{array}$$

where $J_i^{+}$, $J_i^{-}$, and $J_i^{\mathrm{0}}$ are nucleation rates associated with positive, negative, and neutral clusters containing i H₂SO₄ molecules. As a result of scavenging by preexisting particles or wall loss, the steady-state J_i decreases as i increases. To compare with CLOUD measurements, we calculate nucleation at a cluster mobility diameter of 1.7 nm (J_1.7).

Many practical applications require information on the steady state nucleation rates. For each nucleation case presented in this paper, constant values of [H₂SO₄] (i.e., $N_{\mathrm{1}}^{\mathrm{0}}$), [NH₃], T, RH, Q, and $L_i^{+,-,\mathrm{0}}$ are assumed. The preexisting particles with fixed surface area or wall loss serve as a sink for all clusters. Under a given condition, cluster distribution and nucleation rate reach steady state after a certain amount of time. We calculate size-dependent coefficients for a given case and then solve Eqs. (1)–(6) to obtain the steady-state cluster distribution and nucleation rate, with the approach described in Yu (2006b).

Figure 6 Effect of ammonia concentrations ([NH₃]) on effective nucleation rates calculated at a cluster mobility diameter of 1.7 nm (J_1.7, lines) under the stated conditions with two ionization rates (Q) – background ionization, bg (blue), and ionization enhanced by a pion beam, pi (red). Also shown are predictions from the TIMN model, the Atmospheric Cluster Dynamics Code (ACDC) with thermochemistry obtained using the RI-CC2//B3LYP method (McGrath et al., 2012; Kürten et al., 2016), and an empirical parameterization of CLOUD measurements (CLOUDpara) (Dunne et al., 2016) indicated by solid, dashed, and dot-dashed lines, respectively. The symbols refer to CLOUD experimental data (Kirkby et al., 2011; Dunne et al., 2016), with the uncertainties in measured [NH₃] and J_1.7 shown by horizontal and vertical bars, respectively. To be comparable, the CLOUD data points given in Dunne et al. (2016) under the conditions of T = 292 K and RH = 38 % with [H₂SO₄] close to 1.5 × 10⁸ cm⁻³ have been interpolated to the same [H₂SO₄] value (= 1.5 × 10⁸ cm⁻³).

Figure 6 shows a comparison of the model TIMN rates J_1.7 with CLOUD measurements, as a function of [NH₃] under two ionization rates. It should be noted that Dunne et al. (2016) developed a simple empirical parameterization (denoted hereafter as “CLOUDpara”) of binary, ternary, and IIN rates in CLOUD measurements as a function of [NH₃], [H₂SO₄], T, and negative ion concentration. The predictions of CLOUDpara (Dunne et al., 2016) and ACDC based on nucleation thermochemistry obtained using the RI-CC2//B3LYP method (McGrath et al., 2012; Kürten et al., 2016) are also presented in Fig. 6 for comparisons.

Like the CLOUD measurements, the TIMN predictions reveal a complex dependence of J_1.7 on [NH₃], and an analysis of the TIMN results shows this behavior can be explained by the differing responses of negative, positive, and neutral clusters to the presence of ammonia (Fig. 4). Under the conditions specified in Fig. 6, nucleation is dominated by negative ions for [NH₃] < ∼0.5 ppb, by both negative and positive ions for [NH₃] from ∼0.5 ppb to ∼10 ppb (with background ionization) or ∼20 ppb (with pion-enhanced ionization), and by neutrals at higher [NH₃]. According to TIMN, [NH₃] values of at least 0.6–1 ppb are needed before positive ions contribute significantly to nucleation rates – in good agreement with the threshold found in the CLOUD experiments (Kirkby et al., 2011; Schobesberger et al., 2015). TIMN simulations also extend CLOUD data at [NH₃] of ∼1 ppb to include a “zero-sensitivity zone” in the region of 1–10 ppb, followed by a region of strong sensitivity of J_1.7 to [NH₃] commencing at [NH₃] > ∼10–20 ppb. The latter zone may have important implications for NPF in heavily polluted regions, including much of India and China, where [NH₃] may exceed 10–20 ppb (Behera and Sharma, 2010; Meng et al., 2018). It is noteworthy in Fig. 6 that the dependence of J_1.7 on [NH₃] and Q predicted by the ACDC model (McGrath et al., 2012) and the CLOUD data parameterization (Dunne et al., 2016) deviate substantially from the experimental data as well as the TIMN simulations. CLOUDpara does not consider impacts of positive ions and key controlling parameters such as RH and surface area of preexisting particles. Dunne et al. (2016) reported that CLOUDpara is also very sensitive to the approach to parameterize T dependence, showing that the contribution of ternary IIN to NPF below 15 km in altitude has grown from 9.6 % to 37.5 %, after the initial empirical temperature function was replaced with a simpler one.

Figure 7 Comparison of TIMN simulations (solid lines), CLOUDpara predictions (Dunne et al., 2016) (dot-dashed lines), and CLOUD measurements (symbols and data from Dunne et al., 2016) of the dependences of nucleation rates on (a) [H₂SO₄] at five different temperatures (T = 292, 278, 248, 223, and 208 K) and (b) RH at two sets of conditions as specified. [NH₃] is in parts per trillion, by volume. Error bars for the uncertainties in measured [H₂SO₄] (−50 %, +100 %), [NH₃] (−50 %, +100 %), and J_1.7 (overall a factor of 2) are not shown. To be comparable, the CLOUD data points given in Dunne et al. (2016) under the conditions (T, RH, ionization rate) with [NH₃] or [H₂SO₄] close to the corresponding values specified in the figure legends have been interpolated to the same [NH₃] (Fig. 7a) or [H₂SO₄] (Fig. 7b) values.

Figure 7 presents a more detailed comparison of TIMN simulations with CLOUD measurements of J_1.7 as a function of [H₂SO₄], T, and RH. The TIMN model reproduces both the absolute values of J_1.7 and its dependencies on [H₂SO₄], T, and RH in a wide range of temperatures (T = 208–292 K) and [H₂SO₄] (5 × 10⁵–5 × 10⁸ cm⁻³). As expected, nucleation rates are very sensitive to [H₂SO₄] and T. For example, J_1.7 increases by 3 to 5 orders of magnitude with an increase in [H₂SO₄] of a factor of 10 and by roughly 1 order of magnitude for a temperature decrease of 10 ^∘C, except in cases in which the nucleation rate is limited by Q (for example, [H₂SO₄] = ∼10⁸–10⁹ cm⁻³ at T = 278 and 292 K, shown in Fig. 7a). The key difference between CLOUDpara and TIMN predictions is that the dlnJ_1.7 ∕ dln[H₂SO₄] ratio predicted by CLOUDpara is nearly constant while TIMN shows that this ratio depends on both [H₂SO₄] and T. The CLOUD measurements taken at T = 278 K clearly show (in agreement with the TIMN) that dlnJ_1.7 ∕ dln[H₂SO₄] is not constant. CLOUDpara overestimates J_1.7 compared to both measurements and TIMN simulations, except for the case in which T = 278 K and [H₂SO₄] range from ∼7 × 10⁶ to 5 × 10⁷ cm⁻³, with a deviation of CLOUDpara from experimental data and TIMN growing with the lower temperature.

Both CLOUD measurements and TIMN simulations (Fig. 7b) show an important influence of RH on nucleation rates. In particular, CLOUD measurements indicate a 1–5 order of magnitude rise in J_1.7 after RH increases from 10 % to 70 %–80 % and a stronger effect of RH on nucleation rates at higher temperatures under the conditions shown in Fig. 7b. The RH dependence of J_1.7 predicted by the TIMN model is consistent with measurements, being slightly weaker than that measured at high RH.

Figure 8 Model predicted (J_model) versus observed (J_obs) nucleation rates under various conditions of all 377 data points of CLOUD measurements reported in Table S1 of Dunne et al. (2016), with (a) and without (b) the presence of ionization. The data points are grouped according to temperatures as specified in the legend. Vertical error bars show the range of J_model values calculated at 50 % and 200 % of measured [H₂SO₄], corresponding to the uncertainties in measured [H₂SO₄] (−50 %, +100 %). Error bars associated with the uncertainties in measured [NH₃] (−50 %, +100 %) and J_obs (overall a factor of 2) are not shown.

Figure 8 compares TIMN model predictions with all 377 data points of CLOUD measurements reported in data Table S1 of Dunne et al. (2016). The vertical error bars show the range of J_model associated with the uncertainty in the [H₂SO₄] measured (−50 %, +100 %). The effect of uncertainty in measured [NH₃] (−50 %, +100 %) is not included. In the presence of ionization (Fig. 8a), J_model agrees with CLOUD measurements within the uncertainties under mainly all conditions, although J_model tends to be slightly lower than J_obs when T = 292–300 K and J_obs is relatively small (< ∼1 cm⁻³ s⁻¹). For the neutral nucleation (Fig. 8b), the model agrees well with observations at low T (T = 205–223 K) but deviates from observations as T increases. The underprediction of the model for neutral nucleation at T = 278–300 K cannot be explained by the uncertainties in measured [H₂SO₄] and [NH₃]. Apparently for neutral nucleation the model predicts much stronger temperature dependence than the CLOUD measurements. The possible reasons for the difference include the uncertainties in both the model (especially the thermodynamics data and approximation) and measurements. It should be noted that under the conditions of high T and absence of ions, the role of cluster evaporation (i.e., thermodynamics) becomes more important (i.e., higher evaporation and/or generally less tightly bound clusters) and the effect of the possible biases of the used thermochemistry can be more clearly revealed. The contamination (by amines) in the CLOUD measurements (Kirkby et al., 2011) can be another possible reason. The level of contamination in the cloud chamber appears to increase with temperature (Kürten et al., 2016), which may explain the good agreement at low T and increased deviation at higher T. Further research is needed to identify the source of the difference for neutral ternary nucleation at high T.

4 Summary

A comprehensive kinetically based H₂SO₄–H₂O–NH₃ TIMN model, constrained with thermodynamic data from QC calculations and laboratory measurements, has been developed and used to shed new light on physicochemical processes underlying the effect of ammonia on NPF. We show that the stabilizing effect of NH₃ grows with the cluster size and that the reduced effect of ammonia on smaller clusters is caused by weaker bonding that in turn yields lower average NH₃-to-H₂SO₄ ratios. NH₃ was found to impact nucleation barriers for neutral, positively charged, and negatively charged clusters differently due to the large difference in the binding energies of NH₃, H₂O, and H₂SO₄ to small clusters of different charging states. The lowest and highest nucleation barriers are observed in the case of negative ions and neutrals, respectively. Therefore, nucleation of negative ions is favorable, followed by nucleation of positive ions and neutrals. Different responses of negative, positive, and neutral clusters to ammonia result in a complex dependence of ternary nucleation rates on [NH₃]. The TIMN model reproduces both the absolute values of nucleation rates and their dependencies on the key controlling parameters and agrees with the CLOUD measurements for all the cases at the presence of ionization. For the neutral ternary nucleation, the model agrees well with observations at low temperature but deviates from observations as temperature increases.

The TIMN model developed in the present study may be subject to uncertainties associated with the uncertainties in thermodynamic data and interpolation approximation for pre-nucleation clusters. Further measurements and quantum calculations, especially for relatively larger clusters, are needed to reduce the uncertainties. While the TIMN model predicts nucleation rates in a good overall agreement with the CLOUD measurements, its ability to explain the NPF events observed in the real atmosphere is yet to be quantified and will be investigated in further studies.

Appendix A: Quantum-chemical studies of neutral and charged binary and ternary clusters

Thermochemical data for small neutral and charged binary H₂SO₄–H₂O and ternary H₂SO₄–H₂O–NH₃ clusters have been reported in a number of earlier publications (Bandy and Ianni, 1998; Ianni and Bandy, 1999; Torpo et al., 2007; Nadykto et al., 2008; Herb et al., 2011, 2013; Temelso et al., 2012a, b; DePalma et al., 2012; Ortega et al., 2012; Chon et al., 2014; Husar et al., 2014; Henschel et al., 2014, 2016; Kürten et al., 2015). The PW91PW91/6-311$++$G(3df,3pd) method, which is a combination of the Perdue–Wang PW91PW91 density functional with the largest Pople 6-311$++$G(3df,3pd) basis set, has thoroughly been validated and agrees well with existing experimental data. In earlier studies, this method has been applied to a large variety of atmospherically relevant clusters (Nadykto et al., 2006, 2007, 2008, 2014, 2015; Nadykto and Yu, 2007; Torpo et al., 2007; Zhang et al., 2009; Elm et al., 2012, 2013; Leverentz et al., 2013; Xu and Zhang, 2012, 2013; Zhu et al., 2014; Bork et al., 2014; Elm and Mikkelsen, 2014; Peng et al., 2015; Miao et al 2015; Ma et al., 2016) and has been shown to be well suited to study the H₂SO₄–H₂O and H₂SO₄–H₂O–NH₃ clusters, as evidenced by a very good agreement of the computed values with measured cluster geometries, vibrational fundamentals, dipole properties, formation Gibbs free energies (Nadykto et al., 2007, 2008, 2014, 2015; Nadykto and Yu, 2007; Herb et al., 2013; Elm et al., 2012, 2013; Leverentz et al., 2013; Bork et al., 2014), and high-level ab initio results (Temelso et al., 2012a, b; Husar et al., 2012; Bustos et al., 2014).

Table A1 Number of isomers successfully converged at the 6-311 level for selected clusters, along with the enthalpy, entropy, and Gibbs free energy of the most stable isomers (1 hartree = 627.5 kcal mol⁻¹).

Table A2 QC-based stepwise Gibbs free energy change (kcal mol⁻¹) for the addition of one water ($\mathrm{\Delta }{G}_{+\mathrm{W}}^o$), ammonia ($\mathrm{\Delta }{G}_{+\mathrm{A}}^o$), or sulfuric acid ($\mathrm{\Delta }{G}_{+\mathrm{S}}^o$) molecule to form the given positively charged clusters under standard conditions and the corresponding experimental data or semi-experimental estimates.

^a Froyd and Lovejoy (2003a). ^b Meot-Ner (Mautner) et al. (1984). ^c Payzant et al. (1973). ^d Froyd (2002). ^e Froyd and Lovejoy (2012).
^f The $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ values given here were calculated based on experimental $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ values at T = 270 K from Froyd and Lovejoy (2003a) and ΔS values from quantum calculation.

Table A3 Same as Table A2 except for neutral clusters.

^a Nadykto and Yu (2007). ^b Torpo et al. (2007). ^c Ortega et al. (2012). ^d Chon et al. (2007). ^e Kürten et al. (2015).

Table A4 Same as Table A2 except for negatively charged clusters.

^* Froyd and Lovejoy (2003b).

We have extended the earlier QC studies of binary and ternary clusters to larger sizes. The computations have been carried out using a Gaussian 09 suite of programs (Frish et al., 2009). In order to ensure the quality of the conformational search we have carried out a thorough sampling of conformers. We have used both a basin-hopping algorithm, as implemented in Biovia Materials Studio 8.0, and locally developed sampling code. The sampling code is based on the following principle: mesh, with molecules to be added to the cluster placed in the mesh nodes, is created around the cluster, and a blind search algorithm is used to generate the guess geometries. The mesh density and orientation of molecules, as well as the minimum distance between molecules and cluster, are variable. Typically, for each cluster of a given chemical composition a thousand to several thousands of isomers have been sampled. We used a three-step optimization procedure, which includes

• i.

  pre-optimization of initial and guess geometries using the semiempirical PM6 method, separation of the most stable isomers located within 15 kcal mol⁻¹ of the intermediate global minimum and duplicate removal, followed by

• ii.

  optimization of the selected isomers meeting the aforementioned stability criterion by the PW91PW91/CBSB7 method and

• iii.

  the final optimization of the most stable isomers at the PW91PW91/CBSB7 level within 5 kcal mol⁻¹ of the current global minimum using the PW91PW91/6-311$++$G(3df,3pd) method.

Typically, only ∼4 %–30 % of initially sampled isomers reach the second (PW91PW91/CBSB7) level, at which ∼10 %–40 % of isomers optimized with PW91PW91/CBSB7 are selected for the final run. Typically, the number of equilibrium isomers of hydrated clusters is larger than that of unhydrated ones of similar chemical composition. Table A1 shows the numbers of isomers converged at the final PW91PW91/6-311$++$G(3df,3pd) optimization step for selected clusters and HSG (enthalpy, entropy, and Gibbs free energy) values of the most stable isomers used in the present study. The number of isomers optimized at the PW91PW91/6-311$++$G(3df,3pd) level of theory varies from case to case, typically being in the range of ∼10–200.

The computed stepwise enthalpy, entropy, and Gibbs free energies of cluster formation have been thoroughly evaluated and used to calculate the evaporation rates of H₂SO₄ from neutral, positively, and negatively charged clusters.

A1 Positively charged clusters

Table A2 presents the computed stepwise Gibbs free energy changes under standard conditions (ΔG^o) for positive binary and ternary clusters, along with the corresponding experimental data or semi-experimental estimates. Figure 2 in the main text shows ΔG associated with the addition of water ($\mathrm{\Delta }{G}_{+\mathrm{W}}^o$), ammonia ($\mathrm{\Delta }{G}_{+\mathrm{A}}^o$), and sulfuric acid ($\mathrm{\Delta }{G}_{+\mathrm{S}}^o$) to binary and ternary clusters as a function of the cluster hydration number w. Both the absolute values and trends in $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ derived from calculations are in agreement with the laboratory measurements within the uncertainty range of ∼1–2 kcal mol⁻¹ for both QC calculations and measurements. This confirms the efficiency and precision of QC methods in calculating thermodynamic data needed for the development of nucleation models. Nevertheless, it should be noted that the uncertainties in computed free energies of 1–2 kcal mol⁻¹ may lead to large uncertainty in predicted particle formation rates. By increasing or decreasing all Gibbs free energies by 1 kcal mol⁻¹, Kürten at al. (2016) showed that, depending on the conditions, the modeled particle formation rate can change from less than an order of magnitude to several orders of magnitude. Uncertainties estimated by Kürten at al. (2016) represent the upper limit because computed free energies may be overestimated for some clusters and underpredicted for others, which leads to partial or, in some case, full error cancelation.

A2 Neutral clusters

Table A3 presents the computed stepwise Gibbs free energy changes for the formation of ternary S_sA_aW_w clusters under standard conditions. The corresponding binary electrically neutral clusters can be found in previous publications (e.g., Nadykto et al., 2008; Herb et al., 2011). The thermodynamic properties of the S₁A₁ have been reported in a number of computational studies (e.g., Herb et al., 2011; Kurtén et al., 2007; Nadykto and Yu, 2007). However, most of these studies, except for Nadykto and Yu (2007) and Henschel et al. (2014, 2016), did not consider the impact of H₂O on cluster thermodynamics. We have extended the earlier studies of Nadykto and Yu (2007) and Herb et al. (2011) to larger clusters up to S₄A₅ (no hydration) and up to S₂A₂ (hydration included). The free energy of binding of NH₃ to H₂SO₄ (or H₂SO₄ to NH₃) obtained using our method is −7.77 kcal mol⁻¹, which is slightly more negative than values reported by other groups (−6.6–7.61 kcal mol⁻¹) and within less than 0.5 kcal mol⁻¹ of the experimental value of −8.2 kcal mol⁻¹ derived from CLOUD measurements (Kürten et al., 2015).

As it may be seen from Table A3, the NH₃ binding to S₁₋₂W_w weakens as w increases. The average $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ for S₁W_w formation derived from a combination of laboratory measurements and quantum chemical studies are −3.02, −2.37, and −1.40 kcal mol⁻¹ for the first, second, and third hydrations, respectively (Yu, 2007). This indicates that a large fraction of H₂SO₄ monomers in the Earth's atmosphere are likely hydrated. Therefore, the decreasing NH₃ binding strength to hydrated H₂SO₄ monomers implies that RH (and T) will affect the relative abundance of H₂SO₄ monomers containing NH₃. Currently, no experimental data or observations are available to evaluate the impact of hydration (or RH) on $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$. Table A3 shows that the presence of NH₃ in H₂SO₄ clusters suppresses hydration and that $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ for S₂A₂ falls below −2.0 kcal mol⁻¹. This is consistent with earlier studies by our group (Herb et al., 2011) and others (Henschel et al., 2014, 2016) showing that large S_nA_n clusters (n > 2) are not hydrated under typical atmospheric conditions. In the present study, the hydration of neutral S_nA_n clusters at n > 2 is neglected, due to the lack of thermodynamic data.

The number of NH₃ molecules in the cluster (or H₂SO₄ to NH₃ ratio) significantly affects $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ and $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ values. For example, $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ for S₃A_a clusters increases from −7.08 to −16.92 kcal mol⁻¹ and $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ decreases from −16.14 to −8.93 kcal mol⁻¹ as a grows from 1 to 3. For S₄A_a clusters, $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ increases from −7.48 to −16.26 kcal mol⁻¹ and $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ decreases from −17.16 to −11.34 kcal mol⁻¹ as a increases from 2 to 4. $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ for the S₄A₁ cluster is less negative by 1.38 kcal mol⁻¹ than that for S₄A₂. $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ for the S₄A₁ cluster is also quite low (−4.16 kcal mol⁻¹), which might indicate the possible existence of a more stable S₄A₁ isomer, which is yet to be identified. In the presence of NH₃, the uncertainty in the thermochemistry data for S₄A₁ will not significantly affect ternary nucleation rates because most S₄ clusters contain three or four NH₃ molecules.

For the S_sA_a clusters with s=a, $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ increases as the cluster grows while $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ first increases significantly as S₁A₁ is converted into S₂A₂ and then levels off as S₂A₂ is converted into S₄A₄. We also observe a significant drop in $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ when the NH₃ ∕ H₂SO₄ ratio exceeds 1. This finding is consistent with the ACDC model calculation showing that growth of neutral S_sA_a clusters follows the s=a pathway (Schobesberger et al., 2015).

A3 Negative ionic clusters

Table A4 shows ΔG_+W, ΔG_+A, and ΔG_+S needed to form negatively charged clusters under standard conditions, along with available semi-experimental values (Froyd and Lovejoy, 2003b). H₂O binding to negatively charged S⁻S_s clusters significantly strengths with increasing s, from $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ = −0.61 – −1.83 kcal mol⁻¹ at s=1–2 to $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ = −3.5 kcal mol⁻¹ at w=1 and −2.25 kcal mol⁻¹ at w=4 at s=4. $\mathrm{\Delta }{G}_{+\mathrm{W}}^o$ values at s=3 and 4 are slightly more negative (by ∼0.1–0.9 kcal mol⁻¹) than those reported by Froyd and Lovejoy (2003b). Just like H₂O binding, NH₃ binding to S⁻S_s at s < 3 is very weak, with $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ ranging from +2.81 kcal mol⁻¹ at s=0 to −4.85 kcal mol⁻¹ at s=2. However, it significantly increases as s grows. In particular, at $s\ge \mathrm{3}\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ ranges from −11.89 kcal mol⁻¹ for S⁻S₃A₁ to −15.37 kcal mol⁻¹ for S⁻ S₄A₁. NH₃ clearly cannot enter small negative ions. However, it can easily attach to larger negative ions with s≥3, which is consistent with CLOUD measurements (Schobesberger et al., 2015). Since hydration weakens NH₃ binding in S⁻S₃A₁W_w and S⁻S₄A₁W_w clusters, its impacts on the cluster formation and nucleation rates may potentially be important.

In contrast to H₂O and NH₃, binding of H₂SO₄ to small negative ions (s < 3) is very strong. These ions are very stable even when they contain no NH₃ or H₂O molecules. High electron affinity of H₂SO₄ molecules results in the high stability of S⁻S_s at s=1–2. However, the charge effect reduces as s grows. In particular, $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ of S⁻S_s drops from −32.74 kcal mol⁻¹ at s=1 to −10.58 kcal mol⁻¹ and −8.28 kcal mol⁻¹ at s=3 and 4, respectively. At the same time, $\mathrm{\Delta }{G}_{+\mathrm{A}}^o$ increases from 0.08 kcal mol⁻¹ (s=1) to −11.89 kcal mol⁻¹ (s=3) and −15.37 kcal mol⁻¹ (s=4). The hydration of S⁻S_s at s=3 and 4 enhances the strength of H₂SO₄ binding, especially at s=4. $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ values for S⁻S₃₋₄W_w are consistently ∼1.5–3 kcal mol⁻¹ less negative than the corresponding semi-experimental estimates (Table A4). The possible reasons behind the observed systematic difference are yet to be identified and include the use of the low-level ab initio Hartree–Fock method to compute reaction enthalpies and uncertainties in experimental enthalpies in studies by Froyd and Lovejoy (2003b).

NH₃ binding to S⁻S₃ significantly enhances the stability of H₂SO₄ in the cluster by ∼7 kcal mol⁻¹ compared to $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ for the corresponding binary counterpart. The binding of the second NH₃ to S⁻S₃A to form S⁻S₃A₂ is much weaker ($\mathrm{\Delta }{\mathrm{G}}_{+\mathrm{A}}^o$ = −7.27 kcal mol⁻¹) than that of the first NH₃ molecule ($\mathrm{\Delta }{\mathrm{G}}_{+\mathrm{A}}^o=-\mathrm{11.89}$ kcal mol⁻¹). This indicates that most of S⁻S₃A_a can only contain one NH₃ molecule, in a perfect agreement with the laboratory study of Schobesberger et al. (2015). In the case of S⁻S₄, binding of the first ($\mathrm{\Delta }{\mathrm{G}}_{+\mathrm{A}}^{\mathrm{o}}$ = −15.37 kcal mol⁻¹) and second (and −12.23 kcal mol⁻¹) NH₃ molecules to the cluster is quite strong, while the attachment of NH₃ leads to substantial stabilization of H₂SO₄ in the cluster, as evidenced by $\mathrm{\Delta }{G}_{+\mathrm{S}}^o$ growing from −8.28 kcal mol⁻¹ at a=0 to −11.76 kcal mol⁻¹ and −16.71 kcal mol⁻¹ at a=1 and a=2, respectively. The NH₃ binding free energy to S⁻S₄A₂ (to form S⁻S₄A₃) drops to −7.59 kcal mol⁻¹, indicating, in agreement with the CLOUD measurements (Schobesberger et al., 2015), that most S⁻S₄ clusters contain one or two NH₃ molecules.