Impact of non-ideality on reconstructing spatial and temporal variations in aerosol acidity with multiphase buffer theory

Aerosol acidity is a key parameter in atmospheric aqueous chemistry and strongly influences the interactions of air pollutants and the ecosystem. The recently proposed multiphase buffer theory provides a framework to reconstruct long-term trends and spatial variations in aerosol pH based on the effective acid dissociation constant of ammonia ($K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$). However, non-ideality in aerosol droplets is a major challenge limiting its broad applications. Here, we introduced a non-ideality correction factor (c_ni) and investigated its governing factors. We found that besides relative humidity (RH) and temperature, c_ni is mainly determined by the molar fraction of NO${}_{\mathrm{3}}^{-}$ in aqueous-phase anions, due to different NH${}_{\mathrm{4}}^{+}$ activity coefficients between (NH₄)₂SO₄- and NH₄NO₃-dominated aerosols. A parameterization method is thus proposed to estimate c_ni at a given RH, temperature and NO${}_{\mathrm{3}}^{-}$ fraction, and it is validated against long-term observations and global simulations. In the ammonia-buffered regime, with c_ni correction, the buffer theory can reproduce well the $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ predicted by comprehensive thermodynamic models, with a root-mean-square deviation ∼ 0.1 and a correlation coefficient ∼ 1. Note that, while c_ni is needed to predict $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ levels, it is usually not the dominant contributor to its variations, as ∼ 90 % of the temporal or spatial variations in $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ are due to variations in aerosol water and temperature.

1 Introduction

Aerosol acidity strongly influences the thermodynamics and chemical kinetics of atmospheric aerosols and is therefore one essential parameter in evaluating their environmental, health and climate effects (Pye et al., 2020; Zheng et al., 2020). However, direct measurements of aerosol pH in the real atmosphere are not available so far (Pye et al., 2020; Li et al., 2020). The fast equilibrium with ambient air, tiny volume, and high ionic strength and nucleation potential are the main challenges for measurements, especially online or in situ measurements. Several groups are developing new techniques for this purpose (Wei et al., 2018; Craig et al., 2018; Li et al., 2020; Ault, 2020). For example, Wei et al. (2018) developed an in situ Raman microscopy method for pH measurements in microdroplets (diameter ∼ 20 µm), with an uncertainty of ∼ 0.5 pH units. Craig et al. (2018) and Li et al. (2020) developed colorimetric analyses on pH-indicator papers for aerosol pH measurement, which exhibit uncertainties of around 0.4–0.5 pH units. These currently available techniques, however, still need to be developed further for real atmospheric applications.

Due to the lack of direct measurements, modeling tools have been intensively used to calculate the aerosol pH (Fountoukis and Nenes, 2007; Fountoukis et al., 2009; Clegg et al., 2001; Zuend et al., 2008). The results of thermodynamic models are subject to uncertainties in the input parameters (Fountoukis et al., 2009; Pye et al., 2020; Weber et al., 2016; Guo et al., 2016; Guo et al., 2018; Pye et al., 2018; Tao and Murphy, 2019; Hennigan et al., 2015; Guo et al., 2015; Peng et al., 2019). For example, Hennigan et al. (2015) revealed the importance of including gas-phase species in the input, in addition to the full aerosol composition measurements (Fountoukis et al., 2009). Guo et al. (2015) suggests overall uncertainties of ∼ 0.2–0.5 pH units related to aerosol composition. Pye et al. (2020) reviewed major thermodynamic models and showed that the estimated acidity among different models was on average 0.3 pH units but sometimes as much as 1 pH unit.

The recently proposed multiphase buffer theory shows that globally most of the populated urban areas are within the multiphase ammonia-buffered regime (Zheng et al., 2020). In the buffered regions and/or periods, p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ can serve as a proxy of aerosol pH, where $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ is the effective acid dissociation constant of NH₃ in multiphase systems (Sect. 2). Ideally, p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ is fully determined by aerosol water content (AWC) and temperature. However, the non-ideality in aerosols may introduce deviations from the ideal conditions. Here we investigated such deviation and derived a non-ideality correction factor for using p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ as a proxy of aerosol pH. Governing factors of the non-ideality correction factor in aerosol droplets are further explored and discussed, based on which a parameterization method to estimate the non-ideality correction factors is proposed. We also estimated that a constant correction factor of p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ is often good enough to predict pH over a period at a given site or to explain the global pH variations. We thereby provided a way for pH retrieval when chemical measurements are unavailable for the ammonia-buffered regions and periods.

2 Methods

2.1 Effective acid dissociation constant as a proxy of aerosol pH

2.1.1 Acid dissociation constant of NH₃ in bulk solutions, K_a

The definition of acids and bases have been evolving over time. The pioneering Arrhenius theory defined a base as a substance that dissociates in water to form hydroxide (OH⁻) ions (Pfennig, 2015). Therefore, an Arrhenius base can be expressed as B_AOH, which dissociates in water as

$$\begin{array}{}\text{(1a)}& {\mathrm{B}}_{\mathrm{A}}\mathrm{OH}\rightleftharpoons {{\mathrm{B}}_{\mathrm{A}}}^{+}+{\mathrm{OH}}^{-},\end{array}$$

with the corresponding base dissociation constant K_b being

$$\begin{array}{}\text{(1b)}& K_{\mathrm{b}}=\left[{{\mathrm{B}}_{\mathrm{A}}}^{+}\right]\left[{\mathrm{OH}}^{-}\right]/\left[{\mathrm{B}}_{\mathrm{A}}\mathrm{OH}\right].\end{array}$$

In combination with the water dissociation of

$$\begin{array}{}\text{(2)}& {\mathrm{H}}_{\mathrm{2}}\mathrm{O}\rightleftharpoons {\mathrm{H}}^{+}+{\mathrm{OH}}^{-},K_{\mathrm{w}}=\left[{\mathrm{H}}^{+}\right]\left[{\mathrm{OH}}^{-}\right],\end{array}$$

the corresponding acid dissociation constant, defined as $K_{\mathrm{a}}=K_{\mathrm{w}}/K_{\mathrm{b}}$, is thus (Eqs. 2–1)

$$\begin{array}{}\text{(3)}& {\mathrm{B}}_{\mathrm{A}}^{+}+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\rightleftharpoons {\mathrm{B}}_{\mathrm{A}}\mathrm{OH}+{\mathrm{H}}^{+},K_{\mathrm{a}}=\left[{\mathrm{B}}_{\mathrm{A}}\mathrm{OH}\right]\left[{\mathrm{H}}^{+}\right]/\left[{\mathrm{B}}_{\mathrm{A}}^{+}\right].\end{array}$$

The later Brønsted–Lowry theory defined a base as a proton acceptor (Pfennig, 2015) and is expressed as B_BL here. In this sense, an Arrhenius base B_AOH is not regarded as a Brønsted base but rather as salts. The dissociation reaction for a Brønsted base is expressed as

$$\begin{array}{}\text{(4)}& \begin{array}{rl}& {\mathrm{B}}_{\mathrm{BL}}+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\rightleftharpoons {\mathrm{B}}_{\mathrm{BL}}{\mathrm{H}}^{+}+{\mathrm{OH}}^{-},\\ & K_{\mathrm{b}}=\left[{\mathrm{B}}_{\mathrm{BL}}{\mathrm{H}}^{+}\right]\left[{\mathrm{OH}}^{-}\right]/\left[{\mathrm{B}}_{\mathrm{BL}}\right],\end{array}\end{array}$$

and the corresponding K_a is thus (Eqs. 2–4)

$$\begin{array}{}\text{(5)}& {\mathrm{B}}_{\mathrm{BL}}{\mathrm{H}}^{+}\rightleftharpoons {\mathrm{B}}_{\mathrm{BL}}+{\mathrm{H}}^{+},K_{\mathrm{a}}=\left[{\mathrm{B}}_{\mathrm{BL}}\right]\left[{\mathrm{H}}^{+}\right]/\left[{\mathrm{B}}_{\mathrm{BL}}{\mathrm{H}}^{+}\right].\end{array}$$

As NH₃(aq) is actually the water adduct of NH₃, it is often expressed equivalently as NH₃(aq) = NH₃•H₂O(aq) = NH₄OH(aq). In this sense, it can fit in the category of both definitions. In the Arrhenius definition, the base B_AOH = NH₄OH, namely B_A= NH${}_{\mathrm{4}}^{+}$. Therefore, $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}$ is (Eq. 3)

$$\begin{array}{}\text{(6)}& \begin{array}{rl}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}& =\left[{\mathrm{NH}}_{\mathrm{4}}\mathrm{OH}\left(\mathrm{aq}\right)\right]\left[{\mathrm{H}}^{+}\right]/\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\right]\\ & =\left[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)\right]\left[{\mathrm{H}}^{+}\right]/\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\right],\end{array}\end{array}$$

while with the Brønsted–Lowry definition, the base is B_BL=NH₃(aq), and $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}$ is (Eq. 5)

$$\begin{array}{}\text{(7)}& K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}=\left[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)\right]\left[{\mathrm{H}}^{+}\right]/\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\right],\end{array}$$

which is the same as Eq. (6). Therefore, a different definition of bases for the ammonia family (B_A= NH${}_{\mathrm{4}}^{+}$ or B_BL = NH₃(aq)) will lead to the same expression of $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}$, as defined in Zheng et al. (2020). The same applies for other volatile weak bases.

2.1.2 Ideal multiphase acid dissociation constant of NH₃

The multiphase effective acid dissociation constant of NH₃ under ideal conditions, $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}$, depends only on AWC and temperature as (Zheng et al., 2020)

$$\begin{array}{rl}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}& ={\displaystyle \frac{\left[{\mathrm{H}}^{+}\left(\mathrm{aq}\right)\right]\left(\right[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)]+[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)\left]\right)}{\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]}}\\ & =K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}\left(\mathrm{1}+{\displaystyle \frac{{\mathit{\rho }}_{\mathrm{w}}}{H_{{\mathrm{NH}}_{\mathrm{3}}}RT\mathrm{AWC}}}\right),\end{array}$$

where AWC is in micrograms per cubic meters and is mainly determined by air particulate matter concentrations and RH. The ρ_w is water density in micrograms per cubic meter, and AWC $/$ ρ_w represents the aerosol water volume mixing ratio in the air in (m³ water)$/$(m³ air). The [NH₃(g)] represents equivalent molality (in mol kg⁻¹) of gaseous NH₃ in solution (see details in Zheng et al., 2020). The $H_{{\mathrm{NH}}_{\mathrm{3}}}$ is Henry's law constant of NH₃ in mol L⁻¹ atm⁻¹, R is the gas constant of 0.08205 atm L mol⁻¹ K⁻¹, and T is temperature in K.

For typical ambient conditions when AWC varies between 1 to 1000 µg m⁻³, the [NH₃(g)] is usually 10⁵ to 10⁸ times larger than [NH₃(aq)], and the above equation can be simplified into

$$\begin{array}{}\text{(8a)}& \begin{array}{rl}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}& ={\displaystyle \frac{\left[{\mathrm{H}}^{+}\left(\mathrm{aq}\right)\right][{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]}}\\ & =K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}{\displaystyle \frac{{\mathit{\rho }}_{\mathrm{w}}}{H_{{\mathrm{NH}}_{\mathrm{3}}}RT\mathrm{AWC}}}.\end{array}\end{array}$$

Taking negative lognormal on both sides, we have a pH that is related to p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}$ (i.e., −log$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}$) as (Zheng et al., 2020)

$$\begin{array}{}\text{(8b)}& \mathrm{pH}=\mathrm{p}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}+\log {\displaystyle \frac{[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]}}.\end{array}$$

The multiphase buffer capacity of the NH₃–NH${}_{\mathrm{4}}^{+}$ pair reached its local maximum when pH = p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}$, namely when [NH₃(g)] = [NH${}_{\mathrm{4}}^{+}$(aq)]. At a given AWC and T, $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}$ is constant.

2.2 Influences of non-ideality on aerosol pH

For ambient aerosols, the ionic strength (I) is high, and the non-ideality must be considered. Under such non-ideal conditions, the multiphase equilibrium of NH₃ should be expressed as (Zheng et al., 2020)

$$\begin{array}{}\text{(9)}& \begin{array}{rl}& K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}=\\ & {\displaystyle \frac{({\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)}{[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]+{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)}{[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)\left]\right)({\mathit{\gamma }}_{{\mathrm{H}}^{+}}{[\mathrm{H}}^{+}\left(\mathrm{aq}\right)\left]\right)}{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}}\\ & =K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}{\displaystyle \frac{{\mathit{\rho }}_{\mathrm{w}}}{H_{{\mathrm{NH}}_{\mathrm{3}}}RT\mathrm{AWC}}},\end{array}\end{array}$$

where γ_X is the activity coefficient for species X. Note that Eq. (8a) is the simplified expression of Eq. (9) under ideal conditions when all activity coefficients are unity.

Activity coefficients for gases, like ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)}$, are usually treated as unity. Again, for typical ambient conditions [NH₃(g)] is much larger than [NH₃(aq)], and ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)}$[NH₃(aq)] can be omitted. Equation (9) can thus be simplified into

$$\begin{array}{}\text{(10)}& K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}={\displaystyle \frac{{\mathit{\gamma }}_{{\mathrm{H}}^{+}}}{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}}}{\displaystyle \frac{{[\mathrm{H}}^{+}\left(\mathrm{aq}\right){\left]\right[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}}.\end{array}$$

Under non-ideal conditions, pH is usually defined by the proton activity; i.e.,

$$\begin{array}{}\text{(11)}& {\mathrm{pH}}_{\mathrm{a}}=-\log \left({\mathit{\gamma }}_{{\mathrm{H}}^{+}}\right[{\mathrm{H}}^{+}\left]\right).\end{array}$$

However, in thermodynamic models that are most commonly applied in current global models (ISORROPIA II, MOSAIC, etc.), the pH is usually defined as free-H⁺ molality (Pye et al., 2020); i.e.,

$$\begin{array}{}\text{(12)}& \mathrm{pH}={\mathrm{pH}}_F=-\log \left(\right[{\mathrm{H}}^{+}\left]\right).\end{array}$$

The difference in activity- and molality-defined pH (i.e., pH_a and pH_F) has been discussed in a previous study (Pye et al., 2020), which shows that deviations of pH_F from pH_a are larger at lower RH and are usually within 1 unit when RH > 60 % (Pye et al., 2020). To be comparable with results in previous studies, the pH we discuss hereinafter follows the free-H⁺ molality definition. Discussion based on activity-defined pH is detailed in Appendices A and B.

Now we define the multiphase effective acid dissociation constant under non-ideal conditions, $K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$, as

$$\begin{array}{}\text{(13a)}& K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{ni}}={\displaystyle \frac{{[\mathrm{NH}}_{\mathrm{3}}{\left(\mathrm{g}\right)\left]\right[\mathrm{H}}^{+}\left(\mathrm{aq}\right)]}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}},\end{array}$$

which is related to the constant under ideal conditions, $K_{\mathrm{a}}^{\ast ,i}$, by (Eqs. 9, 10, 12 and 13)

$$\begin{array}{}\text{(13b)}& \begin{array}{rl}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{ni}}& ={\displaystyle \frac{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}}{{\mathit{\gamma }}_{{\mathrm{H}}^{+}}}}K{}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}\\ & ={\displaystyle \frac{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}}{{\mathit{\gamma }}_{{\mathrm{H}}^{+}}}}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}{\displaystyle \frac{{\mathit{\rho }}_{\mathrm{w}}}{H_{{\mathrm{NH}}_{\mathrm{3}}}RT\mathrm{AWC}}}.\end{array}\end{array}$$

And the free-H⁺ molality pH is therefore (Eqs. 13a, b)

$$\begin{array}{}\text{(13c)}& \mathrm{pH}=-\log \left(\right[{\mathrm{H}}^{+}\left(\mathrm{aq}\right)])=\mathrm{p}K{}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{ni}}+\log {\displaystyle \frac{[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]}}.\end{array}$$

We now define the non-ideality correction factor c_ni to represent the difference in pH caused by non-ideality. Based on Eqs. (8b) and (13c), c_ni is therefore

$$\begin{array}{}\text{(14a)}& c_{\mathrm{ni}}=\mathrm{p}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{ni}}-\mathrm{p}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,i}.\end{array}$$

And combining Eqs. (13b) and (14a), we have

$$\begin{array}{}\text{(14b)}& c_{\mathrm{ni}}=-\log {\displaystyle \frac{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}}{{\mathit{\gamma }}_{{\mathrm{H}}^{+}}}}.\end{array}$$

Equation (14b) shows the intrinsic determining factors of c_ni, i.e., ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}$ and ${\mathit{\gamma }}_{{\mathrm{H}}^{+}}$. Major influencing factors of c_ni are therefore those influencing the activity coefficients (see Sect. 3.1).

When ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}$ and ${\mathit{\gamma }}_{{\mathrm{H}}^{+}}$ are not available, the c_ni can be alternatively calculated by (Eq. 13a, b)

$$\begin{array}{}\text{(14c)}& \begin{array}{rl}{c}_{\mathrm{ni}}& =\mathrm{p}{K}_{\mathrm{a}}^{\ast ,\mathrm{ni}}-\mathrm{p}{K}_{\mathrm{a}}^{\ast ,\mathrm{i}}\\ & =-\log \left({\displaystyle \frac{{[\mathrm{NH}}_{\mathrm{3}}{\left(\mathrm{g}\right)\left]\right[\mathrm{H}}^{+}\left(\mathrm{aq}\right)]}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}}\right)\\ & +\log \left(K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}{\displaystyle \frac{{\mathit{\rho }}_{\mathrm{w}}}{H_{{\mathrm{NH}}_{\mathrm{3}}}RT\mathrm{AWC}}}\right).\end{array}\end{array}$$

Equation (14c) is valid as [NH₃], [NH${}_{\mathrm{4}}^{+}$] and [H⁺] concentrations will vary as a result of changing c_ni. Note that while [NH₃] $/$ [NH${}_{\mathrm{4}}^{+}$] and pH variations can reflect the c_ni variations and therefore can be used to derive c_ni, they are not the determining factors of c_ni. As shown in Eq. (14b), c_ni is determined by ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}$ and ${\mathit{\gamma }}_{{\mathrm{H}}^{+}}$, which further depends mainly on RH, temperature and the fraction of NO${}_{\mathrm{3}}^{-}$ in anions (see Sect. 3.1).

For some thermodynamics models that predict both the activity coefficients of ions and the gas–particle partitioning of species like the E-AIM model (Sect. 2.3), c_ni can be derived either from Eq. (14b) (activity-based) or Eq. (14c) (gas–particle partitioning based). However, current atmospheric chemical transport models usually adopted the more computation-efficient thermodynamic models (ISORROPIA II, MOSAIC, etc.), in which only the mean activity coefficient of an electrolyte species ij in water, γ_ij, is derived, where i is a cation, while j is an anion (Pye et al., 2020; Fountoukis and Nenes, 2007; Zaveri et al., 2005). For these models, we cannot directly derive ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}$ or ${\mathit{\gamma }}_{{\mathrm{H}}^{+}}$, and c_ni is derived through Eq. (14c) (i.e., from the predicted [NH₃], [NH${}_{\mathrm{4}}^{+}$], [H⁺] and AWC).

2.3 Model simulations

Thermodynamic models. Here we used the E-AIM model (model IV; http://www.aim.env.uea.ac.uk/aim/aim.php, last access: 18 November 2021) (Clegg et al., 1992; Wexler and Clegg, 2002; Friese and Ebel, 2010) to predict both the activity coefficients for individual ions and the gas–particle partitioning. The E-AIM model adopted the Pitzer–Simonson–Clegg model (Clegg et al., 1992, 1998) to calculate single-ion activity coefficients, which includes most comprehensive conditions and has been used as a benchmark (Clegg et al., 1992; Hennigan et al., 2015; Pye et al., 2020). Therefore, both the activity-based pH (pH_a, Eq. 11) and the free-H⁺ molality pH (pH_F, Eq. 12) can be derived (Appendix B). In addition, we also adopted the ISORROPIA v2.3 model (Fountoukis and Nenes, 2007) for comparison, which is computational effective and has been commonly adopted in global and regional models. To reduce the computational cost, the ISORROPIA model calculated only the binary activity coefficients γ_ij using the Kusik–Meissner relationship and Bromley's formula (Fountoukis and Nenes, 2007). Therefore, only the free-H⁺ molality pH (pH_F, Eq. 12) can be derived in ISORROPIA (Appendix D). For example, for a HCl droplet, both γ_H+(aq) and γ_Cl-(aq) are calculated in E-AIM, while only the mean binary activity coefficient of ${\mathit{\gamma }}_{\mathrm{H}\text{-}\mathrm{Cl}}=\sqrt{{\mathit{\gamma }}_{\mathrm{H}+\left(\mathrm{aq}\right)}{\mathit{\gamma }}_{\mathrm{Cl}\text{-}\left(\mathrm{aq}\right)}}$ is estimated in ISORROPIA.

Global models. Spatial variation in c_ni was studied based on the two global models. The global GEOS-Chem model simulations (v11-01) were conducted at a resolution of 2^∘ latitude × 2.5^∘ longitude with 47 vertical layers for 2016. Detailed model settings are provided elsewhere (Zheng et al., 2020). The global EMAC (ECHAM5/MESSy2 for Atmospheric Chemistry) model was run at a resolution of T63 (i.e., ∼ 1.8^∘ × 1.8^∘ at the Equator) with 31 vertical levels for 2016. Detailed EMAC model settings are provided in Appendix C.

3 Results and discussion

3.1 Influencing factors of the non-ideality coefficient

All activity coefficients first depends on RH and temperature. In addition, for ammonium-buffered ambient aerosols, major anions paired with NH${}_{\mathrm{4}}^{+}$ or H⁺ are NO${}_{\mathrm{3}}^{-}$ and SO${}_{\mathrm{4}}^{\mathrm{2}-}$. The ratio of mean activity coefficients is therefore expected to differ when they are mainly combined with SO${}_{\mathrm{4}}^{\mathrm{2}-}$ (i.e., ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}{\mathrm{HSO}}_{\mathrm{4}}}/{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}}$) or NO${}_{\mathrm{3}}^{-}$ (i.e., ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}{\mathrm{NO}}_{\mathrm{3}}}/{\mathit{\gamma }}_{{\mathrm{HNO}}_{\mathrm{3}}}$).

Figure 1 shows the dependence of c_ni under different systems (Appendix A), as predicted by the gas–particle portioning (Eq. 14c) with E-AIM (Fig. 1a, c, e) and ISORROPIA II (Fig. 1b, d, f). Based on both models, c_ni differs much between the NH₃–H₂SO₄ system (Fig. 1a, b) and the NH₃–HNO₃–H₂SO₄ system (Fig. 1c, d), even at the same RH and temperature. The difference is still large when compared at the same ionic strength and temperature (Fig. A1), illustrating that the difference is mainly due to the ion-pair specific binary activity coefficients, ${\mathit{\gamma }}_{ij}^o$ (Zaveri et al., 2005; Fountoukis and Nenes, 2007; Clegg et al., 1992) (Appendix B; Fig. B1).

Figure 1 The non-ideality correction factor, c_ni, estimated by E-AIM (a, c, e) and ISORROPIA (b, d, f) for different aerosol systems. (a, b) NH₃–H₂SO₄–H₂O system with aerosols dominated by (NH₄)₂SO₄ at varying RH and temperature conditions; (c, d) NH₃–HNO₃–H₂O system with aerosols dominated by NH₄NO₃ at varying RH and temperature conditions, and (e, f) Na⁺–NH₃–HNO₃–H₂SO₄–H₂O system with varying chemical profiles at 288.15 K and RH of 73 %. The chemical profiles in (e) and (f) are characterized by the fraction of NO${}_{\mathrm{3}}^{-}$ in anions(aq) and NVCs $/$ anions(aq), where the non-volatile cations (NVCs) are assumed to be Na⁺ only here. The assumed RH and T conditions in (e) and (f) are marked as blacked stars in (a)–(d), while the chemical profiles of (a)–(d) are marked by the corresponding letter in (e) and (f), respectively.

Due to the large difference in c_ni between NH₄NO₃ and (NH₄)₂SO₄-dominated aerosols, the c_ni at a given RH and temperature conditions is therefore sensitive to the anion profiles, as characterized by the fraction of NO${}_{\mathrm{3}}^{-}$ in anions(aq), $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$, of

$$\begin{array}{}\text{(15a)}& f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}\left(\mathrm{\mu }\mathrm{eq}/\mathrm{\mu }\mathrm{eq}\right)=\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)\right]/\left[\mathrm{anions}\left(\mathrm{aq}\right)\right],\end{array}$$

$$\begin{array}{}\text{(15b)}& \left[\text{anions(aq)}\right]=\mathrm{2}\left[{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}\left(\mathrm{aq}\right)\right]+\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)\right]+\left[{\mathrm{Cl}}^{-}\left(\mathrm{aq}\right)\right].\end{array}$$

The $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ is proportional to NO${}_{\mathrm{3}}^{-}/$SO${}_{\mathrm{4}}^{\mathrm{2}-}$ molar ratios when Cl⁻ is negligible. In comparison, the cation profiles, or the relative abundances of non-volatile cations (NVCs; total cations from Na⁺, Ca²⁺, K⁺ and Mg²⁺), play a minor role as their influence is more indirect (Fig. 1e, f).

3.2 Comparison of c_ni estimated by E-AIM and ISORROPIA

As discussed in Sect. 2.2, for E-AIM c_ni can be estimated either by activity coefficients (Eq. 14b) or gas–particle portioning (Eq. 14c), and the results agreed perfectly (black lines in Fig. 2). Therefore, the c_ni estimation with E-AIM is calculated by the gas–particle portioning (Eq. 14c) hereinafter, the same as for ISORROPIA.

Figure 2 Comparison of different c_ni estimation methods for three representative aerosol systems. The c_ni are compared under the same conditions (i.e., same RH, temperature and chemical profiles). The x values are c_ni estimated by activity coefficients (Eq. 14b) with the E-AIM model, and the y values include c_ni estimated by gas $/$ particle ratios (Eq. 14c) with the E-AIM (black lines) and ISORROPIA (blue dots) models. The systems are the same as in Fig. 1.

Although showing the same influencing factors, c_ni estimated by E-AIM and ISORROPIA is not identical (Fig. 1). Especially for the NH₃–H₂SO₄–H₂O system (i.e., (NH₄)₂SO₄-dominated aerosols), E-AIM (Fig. 1a) and ISORROPIA (Fig. 1b) even predicted reversed trends in c_ni dependence on RH and temperature. This is more clearly shown in Fig. 2 (blue dots), where c_ni by E-AIM and ISORROPIA in the same conditions (i.e., same RH, temperature and chemical profiles) is compared. As shown in Fig. 2a, while c_ni predicted by E-AIM ranged −0.3 to 0.5 for (NH₄)₂SO₄-dominated aerosols, that predicted by ISORROPIA was always larger than 0.1. This is mainly caused by the difference in calculated activity coefficients between ISORROPIA and E-AIM (Eq. 14b; see details in Appendix D, Figs. D1 and D2).

Despite the large difference in predicted c_ni for the NH₃–H₂SO₄–H₂O system, the E-AIM and ISORROPIA models generate similar predictions of AWC and therefore a similar ideal constant of $K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}$ (Fig. D1a). Combined with different c_ni, this would lead to a different prediction of [H⁺(aq)][NH₃(g)]$/$[NH${}_{\mathrm{4}}^{+}$(aq)] by the two models (Eq. 14c). However, with the constraint of charge balance and mass conservations of ammonia (Appendix D), the disagreement in the predicted molar ratios of NH₃(g)$/$NH${}_{\mathrm{4}}^{+}$(aq) between these two models is relatively small (4 %–6 %; Fig. D1b), and most of the c_ni variations are allocated to the [H⁺], or pH, predictions (Fig. D1c).

Unlike the NH₃–H₂SO₄–H₂O system, c_ni estimated by ISORROPIA generally agrees well with E-AIM when HNO₃ is present in the system (Fig. 2b, c), although it tends to be somewhat higher than E-AIM. This indicates that constraints from NH₃–HNO₃ equilibriums are quite important in estimating c_ni with ISORROPIA (see details in Appendix D). Under ambient conditions, there is barely places with negligible HNO₃, thus the ISORROPIA-predicted c_ni generally agreed with E-AIM (Sect. 3.4).

With the known governing factors, here we propose a parameterization method to estimate c_ni at a given RH, temperature and $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$, with lookup tables generated by comprehensive thermodynamic models, E-AIM and ISORROPIA (“AIM_molality” database and “ISORROPIA_molality” database as in Data S1). In addition, the parameterized c_ni for activity-based pH (Eq. 11; Appendix B; Fig. B1) is also available (“AIM_activity” database in Data S1). A Matlab code to get c_ni is also provided (Data S1). Example slices of this c_ni parameterization based on AIM_molality estimations are shown in Fig. 3. Note that this parameterization method is aimed at the NH₃–HNO₃–H₂SO₄–H₂O system, assuming no NVCs. We will show that this assumption is acceptable under most cases in the following sections.

Figure 3 Example slices of the c_ni parameterization based on “AIM_molality” estimations as given in Data S1 in the Supplement.

3.3 Validation and applications with long-term observations

To validate the c_ni parameterization method under actual ambient conditions, we here show an example application based on the long-term measurements in Toronto (Tao and Murphy, 2019) (Fig. 4). From 2007 to 2016, Toronto resides in the ammonia-buffered regime for ∼ 80 % of the time, and the model-predicted pH based on the measured chemical compositions follows nicely with the variation in actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ estimated by thermodynamic models (Eq. 14), for both E-AIM (Fig. 4a) and ISORROPIA (Fig. 4c). Parameterized c_ni agreed quite well with the actual ones for both models (Fig. 4b, d, black circles), with R² being 0.99 in both cases and the both corresponding root-mean-square deviations (RMSDs) being ∼ 0.1.

Figure 4 Comparison of pH and actual and estimated p$K_{\mathrm{a},\mathrm{ni}}^{\ast }$ based on the 10-year observations in Toronto. Data were taken from Canada's National Air Pollution Surveillance (NAPS) Program, as detailed in Tao and Murphy (2019). Predications are based on (a–b) the E-AIM model and (c–d) the ISORROPIA model. The “parameterized” series in (b, d) are predicted by the parameterization method proposed with the input of the observed RH, temperature and model-predicted fraction of nitrates in anions. The annual mean and monthly mean are based on mean c_ni of an arbitrary example year of 2012.

Figure 4 also suggests that most of the variation in actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ comes from the variation in ideal constants (p$K_{\mathrm{a}}^{\ast ,\mathrm{i}}$), not the non-ideality. For example, now we assume the full aerosol and gas measurements were conducted only in a calibration year of 2012, based on which the annual mean and monthly mean c_ni can be derived (Fig. E1). Annual mean c_ni is 0.4 for E-AIM and 0.8 for ISORROPIA estimations. When we use the annual mean c_ni as a constant correction (i.e., estimated p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}=$ p$K_{\mathrm{a}}^{\ast ,\mathrm{i}}+$ annual mean), fluctuation in the estimated p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ would actually all come from p$K_{\mathrm{a}}^{\ast ,\mathrm{i}}$. However, this estimated p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ can already explain ∼ 90 % of the variations in actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ (red dots in Fig. 4b, d), illustrating the dominance of p$K_{\mathrm{a}}^{\ast ,\mathrm{i}}$ (i.e., AWC and temperature fluctuations) over non-ideality. In comparison, applying the month-dependent c_ni values (blue dots in Fig. 4b, d) makes little difference with the annual constant estimations (R² differed only by 1 %).

Figures 4 and E1 illustrate that a constant c_ni is often good enough at a given site. Full aerosol species measurements for a whole year, or under periods representative of annual-average conditions (like spring or fall seasons for Toronto; Fig. E1), are recommended in determining the localized c_ni, which, together with AWC and temperature measurements, could already provide a good approximation of the aerosol pH. This is especially useful in retrieving the acidity variations when full chemical measurements are not available in the long run.

3.4 Validation and application against global model simulations

We further investigated the influence of non-ideality in explaining the spatial variations in aerosol acidity based on global model simulations. On the global scale, the fraction of NO${}_{\mathrm{3}}^{-}$ in aqueous-phase anions depends on two factors: the total nitrate (gas + particle phase) to sulfate ratios and the partitioning of total nitrates. When total nitrate ≪ sulfates, the aerosols would be dominated by (NH₄)₂SO₄ even if all the nitrates are partitioned into the particle phase. In this case, a non-ideality correction factor can be estimated from Fig. 1a and b at known RH and temperature. However, both GEOS-Chem and the EMAC results show that this criterion is barely met for the ammonia-buffered regions. Besides, for all the reported observation results we know of, only the summertime southeastern US (Weber et al., 2016) has a total nitrate that is < 5 % of the sulfate (charge ratios). Therefore, under most conditions, c_ni largely depends on the partitioning of total nitrates, and an estimation of $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ is needed to derive the correction factor.

Figure 5 Comparison of actual and estimated p$K_{\mathrm{a},\mathrm{ni}}^{\ast }$ based on the GEOS-Chem global simulations in 2016. Predications are based on (a–b) the E-AIM model and (c–d) the ISORROPIA model. The “parameterized” series are based on the parameterization method proposed in this study, while the global means are based on mean c_ni calculated from thermodynamic models under each scenario.

Figure 5 shows the estimated p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ against actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ based on GEOS-Chem simulations, and that based on EMAC simulations is shown in Fig. C1. Three scenarios are assumed to examine the sensitivity of p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ prediction with c_ni values: (a) constant temperature (T) of 288 K and RH of 73 %, (b) constant RH of 73 % but with annual-average temperatures for each site, and (c) annual-average T and RH for each site. For all scenarios, annual mean chemical compositions for the ammonia-buffered surface regions (Zheng et al., 2020) are used, and c_ni is estimated by both E-AIM and ISORROPIA II models. Similar to Fig. 4, in the “parameterized” series c_ni is estimated by the parameterization method proposed in this study with RH, T and $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ at certain model grids, while in the “global mean” series, c_ni is assumed to be constant as the average of actual c_ni estimated by the thermodynamic models under each scenario, which is ∼ 0.6 for the E-AIM model and ∼ 0.8 for the ISORROPIA model.

Based on GEOS-Chem simulations, the parameterized c_ni (black dots in Fig. 5) work nicely in reproducing actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$, with R² near 1 under all scenarios and the RMSD of < 0.03 for the ISORROPIA model and ∼ 0.1 for the AIM model. Again, we found that variations in c_ni are much smaller than the variation in p$K_{\mathrm{a}}^{\ast ,\mathrm{i}}$ caused by particulate matter concentrations and temperatures. With a constant global-mean c_ni correction (i.e., assuming a global average $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$) (blue dots in Fig. 5), the estimated p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ can already explain over 93 % of the variations in actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$, with/without considering the influence of meteorology on non-ideality alike. Correspondingly, it can already explain ∼ 70 % of the aerosol pH variations (Zheng et al., 2020), where the pH is further subject to variations in NH₃(g) and NVCs (Eq. 8; Zheng et al., 2020).

The EMAC simulations show similar patterns with GEOS-Chem results. Estimated p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ with the parameterized c_ni corrections agreed well with actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$, with R² over 0.94 for the E-AIM model and over 0.91 for the ISORROPIA model (Fig. C1). This is somewhat lower than the Toronto site (Fig. 4) or the GEOS-Chem result (Fig. 5), which is due to the larger variations in the simulated chemical profiles (e.g., importance of NVCs and Cl⁻). The constant c_ni assumption (blue dots in Fig. C1) works similarly with the parameterized ones when the influence of meteorology is excluded (Fig. C1a, d) or when spatial variations in temperatures are considered (Fig. C1b, e). When spatial variations in both temperature and RH (Fig. C1c, f) are considered, the constant c_ni assumption works worse than the parameterized ones but is still acceptable (R² being 0.75 for E-AIM and 0.69 for ISORROPIA).

Note that under all conditions, the global mean method tends to overestimate c_ni when the actual p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ of NH₃ is smaller than 2 (Figs. 5, C1). That is caused by $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$. The low p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ indicates low AWC levels (Zheng et al., 2020) and relatively low pH levels (Eq. 13). Under such conditions, HNO₃ tends to stay in the gas phase (Nenes et al., 2020), corresponding to a low $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ of ∼ 0. In comparison, the global-mean c_ni corresponds to the global-mean simulated $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ of ∼ 0.4. As c_ni increases with increasing $f_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ (Fig. 1e, f), the global-mean c_ni would tend to overestimate actual low p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ conditions (i.e., < ∼ 2).

4 Conclusions

Overall, we found that the non-ideality correction is needed for using p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ of NH₃ as a proxy of aerosol pH in ammonia-buffered regimes. This correction factor, c_ni, generally ranges from 0.3–1.1 and mainly depends on RH, temperature and the fraction of nitration in aqueous-phase anions. E-AIM generally predicted a lower c_ni than the ISORROPIA model. We proposed a parameterization method to estimate the c_ni, which works quite well, as validated against both long-term observations and global simulations. Although the correction is needed in estimating the ammonia p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ levels, the variations in p$K_{\mathrm{a}}^{\ast ,\mathrm{ni}}$ is often much less sensitive to the non-ideality than to aerosol water content and temperature. Therefore, a constant correction factor of p$K_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast }$ is often good enough to predict pH over a period at a given site, or to explain the global pH variations. We thereby provided a way for pH retrieval when chemical measurements are unavailable for the ammonia-buffered regions.

Appendix A: Scenario settings for different systems

In Figs. 1, 2, A1, B1, D1 and D2, we assumed three systems, with the settings detailed below.

Figure A1 Ionic strength (I) and the non-ideality correction factor, c_ni, as calculated by E-AIM (a, c, e) and ISORROPIA (b, d, f) under different aerosol systems. The systems are the same as in Fig. 1, while the RH in Fig. 1a–d and c_ni in Fig. 1e–f are replaced by I.

NH₃–H₂SO₄–H₂O system. For this system, we assumed a constant input with 0.5 µmol m⁻³ of total sulfate (i.e., 1 µeq m⁻³of anions) and 2 µmol m⁻³ of total ammonia. This ratio is to ensure that the system pH is around the maximum buffering capacity of ammonia. However, we found that for the ISORROPIA model, the solver with only ammonia and sulfate inputs is not stable, with the predicted pH often larger than 7. We thereby introduced 0.015 µmol m⁻³ of Na⁺ (3 % of the total sulfate molar concentrations or 1.5 % of the anions), which exerted little influence on the ionic environments (difference in E-AIM results less than 3 %) but will change the ISORROPIA subroutine solver called. The RH and temperature are then varied at different values to check the influence.

NH₃–HNO₃–H₂O system. For this system, we assumed a constant input with 1 µmol m⁻³ of total nitrate (also 1 µeq m⁻³of anions) and 2 µmol m⁻³ of total ammonia and then varied the RH and temperatures to derive non-ideality correction factors.

Na⁺–NH₃–HNO₃–H₂SO₄–H₂O system. For this system, we fixed the RH at 73 % and temperature at 288.15 K, 2 µmol m⁻³ of total ammonia and a fixed concentration of total anions of 1 µeq m⁻³. The nitrate $/$ sulfate ratios are then varied (but keeping their total charges the same) to get different nitrate fractions. For example, when the input sulfate is 0.25 µmol m⁻³ equalling 0.5 µeq m⁻³ of anions, the input total nitrate is then set to 0.5 µmol m⁻³, corresponding to a total anion of the system of 1 µeq m⁻³. Meanwhile, the ratio of NVCs (here assumed to be Na⁺ only) to anions is also varied and combined with the different nitrate $/$ sulfate ratios to generate different simulation conditions.

Appendix B: Non-ideality correction factor for activity-based pH definitions

Figure B1 Dependence of the non-ideality correction factor for activity-based pH definitions, c_nia (i.e., −log(${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}}$), as estimated by E-AIM. (a) NH₃–HNO₃–H₂O system with aerosols dominated by NH₄NO₃ at varying RH and temperature conditions and (b) Na⁺–NH₃–HNO₃–H₂SO₄–H₂O system with varying chemical profiles at 288.15 K and a RH of 73 %. Note that the c_nia for the NH₃–H₂SO₄–H₂O system (i.e., (NH₄)₂SO₄-dominated aerosols) is not shown, as it varied little (ranging 0.44–0.47) over all the RH and temperature ranges explored.

With an activity-based pH definition (i.e., pH $=-\log$(${\mathit{\gamma }}_{{\mathrm{H}}^{+}}$ [H⁺]), the multiphase buffer theory can be rewritten as

$$\begin{array}{}\text{(B1a)}& {\displaystyle }\begin{array}{rl}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{nia}}& ={\displaystyle \frac{{[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}}({\mathit{\gamma }}_{{\mathrm{H}}^{+}}{[\mathrm{H}}^{+}\left(\mathrm{aq}\right)\left]\right)={\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}K{}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{i}}\\ & ={\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}\left(\mathrm{1}+{\displaystyle \frac{{\mathit{\rho }}_{\mathrm{w}}}{H_{{\mathrm{NH}}_{\mathrm{3}}}RT\mathrm{AWC}}}\right),\end{array}\text{(B1b)}& {\displaystyle }\begin{array}{rl}\mathrm{pH}& =-\log \left({\mathit{\gamma }}_{{\mathrm{H}}^{+}}{[\mathrm{H}}^{+}\left(\mathrm{aq}\right)]\right)\\ & =\mathrm{p}{K}_{\mathrm{a},{\mathrm{NH}}_{\mathrm{3}}}^{\ast ,\mathrm{nia}}+\log {\displaystyle \frac{[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]}},\end{array}\end{array}$$

where $K_{\mathrm{a}}^{\ast ,\mathrm{nia}}$ is the multiphase effective acid dissociation constant under non-ideal conditions. The difference in pH caused by non-ideality, c_nia, is therefore

$$\begin{array}{}\text{(B2)}& c_{\mathrm{nia}}=\mathrm{p}{K}_{\mathrm{a}}^{\ast ,\mathrm{nia}}-\mathrm{p}{K}_{\mathrm{a}}^{\ast ,\mathrm{i}}=-\log {\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}.\end{array}$$

That is, the non-ideality correction factor for activity-based pH is actually the ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}$, which can be calculated with the more comprehensive models like E-AIM. The E-AIM calculated mole-fraction-based activity coefficient (f_i) that can be converted to the molality-based activity coefficient (γ_i) by (Pye et al., 2020; Peng et al., 2019)

$$\begin{array}{}\text{(B3)}& {\mathit{\gamma }}_i=f_i{x}_{\mathrm{w}}=f_i{x}_i/\left(m_i{M}_{\mathrm{w}}\right),\end{array}$$

where x_i and m_i are, respectively, the mole fraction and molality of species i and x_w and M_w are, respectively, the mole fraction and molecular weight of water. All these variables are given in E-AIM outputs. Major influencing factors of c_nia are also RH, temperature and the fraction of NO${}_{\mathrm{3}}^{-}$ in anions in the aqueous phase (aq), as shown in Fig. B1.

Appendix C: EMAC model settings

In this section, we will only focus on the model settings for EMAC simulations, while for the GEOS-chem model settings, please refer to Zheng et al. (2020). We used the global ECHAM/MESSy Atmospheric chemistry – Climate (EMAC) model, which is a numerical chemistry and climate simulation system that includes submodels describing tropospheric and middle-atmosphere processes and their interaction with oceans, land and human influences (Jöckel et al., 2010). The core atmospheric model is the 5th generation European Centre Hamburg general circulation model (ECHAM5) (Roeckner et al., 2006), which has been modularized and which improved submodels and updates of boundary layer, radiation, cloud and convection routines have been introduced. The EMAC model development is coordinated within an international consortium: see https://www.messy-interface.org (last access: 18 November 2021). For the present study we applied EMAC (ECHAM5 version 5.3.02, MESSy version 2.54.0) in the T63L31 resolution, i.e., with a spherical truncation of T63 (corresponding to a quadratic Gaussian grid of approximately 1.8 by 1.8^∘ in latitude and longitude) with 31 vertical hybrid terrain-following pressure levels up to 10 hPa in the lower stratosphere. Meteorological conditions as in ERA-Interim data from the European Centre for Medium- range Weather Forecasts (ECMWF) were simulated by the model by applying a “nudging” technique (Jöckel et al., 2006). EMAC simulates gas-phase and heterogeneous chemistry through the MECCA submodel, which accounts for the photochemical oxidation of natural and anthropogenic emissions, including a comprehensive account of volatile organic carbon compounds (Sander et al., 2019). Aerosol microphysical processes and gas–particle partitioning are simulated with the GMXe submodel (Pringle et al., 2010; Pozzer et al., 2012), which describes the aerosol size distribution by seven interacting lognormal modes (four hydrophilic and three hydrophobic). The aerosol composition can vary between these modes (externally mixed) and is uniform within each mode (internally mixed). The inorganic aerosol composition is computed with the ISORROPIA-II thermodynamic equilibrium submodel (Fountoukis and Nenes, 2007). It calculates the gas–liquid–solid equilibrium partitioning of inorganic compounds and water. The composition and atmospheric evolution of organic aerosol compounds is simulated with the ORACLE submodel, which represents volatility classes of organics through their effective saturation concentrations (Tsimpidi et al., 2018). For this work the anthropogenic emissions from EDGAR (Emissions Database for Global Atmospheric Research v4.3.2) (Crippa et al., 2018) were applied, as well as the GFAS (Global Fire Assimilation System, v1.0) (Kaiser et al., 2012) for biomass burning emissions. The EMAC results are shown in Fig. C1.

Figure C1 Same as Fig. 5, but based on EMAC results.

Appendix D: Potential reasons for discrepancies in predicting aerosol pH by ISORROPIA and E-AIM for the NH₃–H₂SO₄–H₂O system

Figure D1 Drivers of the difference in c_ni estimated by the ISORROPIA and E-AIM models for the NH₃–H₂SO₄–H₂O system. The ${\mathit{\zeta }}_{\mathrm{g}/\mathrm{p}}$ of NH₃ indicates the molar ratios of NH₃(g) to particle-phase NH${}_{\mathrm{4}}^{+}$.

In this study, we applied ISORROPIA version 2.1 (Fountoukis and Nenes, 2007) and E-AIM (model IV; http://www.aim.env.uea.ac.uk/aim/aim.php, last access: 18 November 2021) (Clegg et al., 1992; Wexler and Clegg, 2002; Friese and Ebel, 2010), and the following description and discussion refer to these versions of the models. For the NH₃–H₂SO₄–H₂O system, we found that by assuming the same input of 0.5 µmol m⁻³ of total sulfate (i.e., 1 µeq m⁻³of anions) and 2 µmol m⁻³ of total ammonia and varying the RH (60 %–90 %) and temperatures (265–300 K), ISORROPIA predicted a very high aerosol pH of about 13 (12.6–13.2), while the E-AIM-predicted pH ranged from 2–5, which is obviously more realistic. However, by introducing only an small amount of Na⁺ (0.015 µmol m⁻³ or 3 % of the total sulfate), the ISORROPIA-predicted pH dropped dramatically to 2–5 (Fig. D1), while the E-AIM-predicted pH changed little compared to the no-Na⁺ predictions (pH increased systematically by 0.03 with both R² and slope being 1). Besides, the predicted pH assuming only HNO₃ and NH₃ inputs (NH₃–HNO₃–H₂O system) agreed well between ISORROPIA and E-AIM.

Figure D2 Comparison of activity coefficients for different species. (a–c) Comparison of activity coefficients involved in ISORROPIA A2 subcase calculations, as predicted by ISORROPIA and E-AIM. (d) Comparison of activity coefficients involved in ISORROPIA D3 subcase calculations, as predicted by ISORROPIA and E-AIM. (e) Mean activity coefficients predicted by ISORROPIA that are involved in ISORROPIA G5 subcase calculations.

We found that the dramatic changes in ISORROPIA-predicted pH levels with or without small amounts of Na⁺ and NO${}_{\mathrm{3}}^{-}$ additions are related to the different calculation procedures among subcases. Here we focused on subcases under the metastable and sulfate-poor (i.e., total potential cations, including total ammonia ([NH₃]_t) and NVCs, exceed twice the molar ratios of total sulfate ([H₂SO₄]_t)) conditions.

In ISORROPIA, when there are only NH₃ and H₂SO₄ (i.e., the “pure” NH₃–H₂SO₄–H₂O system), the corresponding subcase is “A2”. As detailed below, for this subcase, activity coefficients included in the final calculations are ${\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}}$, ${\mathit{\gamma }}_{\mathrm{2}\mathrm{H}-{\mathrm{SO}}_{\mathrm{4}}}$ and ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}-{\mathrm{HSO}}_{\mathrm{4}}}$. As shown in Fig. D2a–c, for all these three values, there is a large difference between E-AIM and ISORROPIA estimations (note that log scales are used for ${\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}}$ and ${\mathit{\gamma }}_{\mathrm{2}\mathrm{H}-{\mathrm{SO}}_{\mathrm{4}}}$ plots). Therefore, it is not surprising that there is a large discrepancy between the predicted pH from subcase A2 of ISORROPIA and E-AIM.

In comparison, the subcase would change to “D3” when HNO₃ is introduced to the system. As detailed below, for this subcase, only ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}}$ is involved in the calculations. As shown in Fig. D2d, although ISORROPIA still shows a different trend than E-AIM, it is, however, at least on the same order of magnitude as the one predicted by E-AIM.

By introducing a small amount of Na⁺ into the NH₃–H₂SO₄–H₂O system, the calculation procedure of ISORROPIA would change from A2 to G5 (an Na⁺–NH₃–H₂SO₄–HNO₃–HCl–H₂O aerosol system). For the G5 subcase, we noticed two issues: (1) although the total HNO₃ is zero, the model still tried to predict ${\mathit{\gamma }}_{\mathrm{H}-{\mathrm{NO}}_{\mathrm{3}}}$ and ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}}$; (2) as it was using Cl⁻ as the x variable for the final solutions, a small amount of Cl⁻ was always present, which was introduced by the model so the calculation procedures could go on. The relevant values are shown in Fig. D2e. In comparison, E-AIM predicted no NO${}_{\mathrm{3}}^{-}$ or Cl⁻, and the activity coefficients of other relevant species change little with the no-Na⁺ case. Therefore, we could not perform a comparison between ISORROPIA and E-AIM for this case (as there is no ${\mathit{\gamma }}_{{\mathrm{NO}}_{\mathrm{3}}}$ or γ_Cl in E-AIM). Based on the pH and non-ideality comparisons (Fig. D1), however, we could see that the NH₃ partitioning estimated in this way is far more realistic than the A2 subcase.

Calculation principles for subcase A2 (an NH₃–H₂SO₄–H₂O aerosol system). For the subcase A2, the major constraining equations include the [SO${}_{\mathrm{4}}^{\mathrm{2}-}$]$/$[HSO${}_{\mathrm{4}}^{-}$] equilibriums, gas–particle partitioning of ammonia and charge balance:

$$\begin{array}{}\text{(D1a)}& {\displaystyle }\begin{array}{rl}& {\mathrm{HSO}}_{\mathrm{4}}^{-}\rightleftharpoons {\mathrm{H}}^{+}+{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-},\\ & K_{\mathrm{a},{\mathrm{HSO}}_{\mathrm{4}}}={\displaystyle \frac{[{\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}\left(\mathrm{aq}\right)]}{{[\mathrm{HSO}}_{\mathrm{4}}^{-}\left(\mathrm{aq}\right)]}}{\displaystyle \frac{{\mathit{\gamma }}_{{\mathrm{H}}^{+}\left(\mathrm{aq}\right)}{\mathit{\gamma }}_{{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}\left(\mathrm{aq}\right)}}{{\mathit{\gamma }}_{{\mathrm{HSO}}_{\mathrm{4}}^{-}\left(\mathrm{aq}\right)}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}\left(\mathrm{aq}\right)]}{{[\mathrm{HSO}}_{\mathrm{4}}^{-}\left(\mathrm{aq}\right)]}}{\displaystyle \frac{{\mathit{\gamma }}_{\mathrm{2}\mathrm{H}\text{-}{\mathrm{SO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{3}}}{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}},\end{array}\text{(D1b)}& {\displaystyle }\begin{array}{rl}& {\mathrm{NH}}_{\mathrm{4}}^{+}\rightleftharpoons {\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)+{\mathrm{H}}^{+},\\ & K_{\mathrm{ag},\mathrm{NH}\mathrm{3}}=K_{\mathrm{w}}/K_{\mathrm{bg},{\mathrm{NH}}_{\mathrm{3}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]{\mathit{\gamma }}_{{\mathrm{H}}^{+}\left(\mathrm{aq}\right)}}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}}{\displaystyle \frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}},\end{array}\text{(D1c)}& {\displaystyle }\left(\right[{\mathrm{NH}}_{\mathrm{4}}^{+}]+[{\mathrm{H}}^{+}\left]\right)/\left(\mathrm{2}\right[{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}]-[{\mathrm{HSO}}_{\mathrm{4}}^{-}\left]\right)-\mathrm{1}=\mathrm{0}.\end{array}$$

With these three equations and known total ammonia ([NH₃]_t) and total sulfate ([H₂SO₄]_t), we have

$$\begin{array}{}\text{(D2)}& \begin{array}{rl}& \left[{\mathrm{NH}}_{\mathrm{3}}{]}_{\mathrm{t}}{C}_{\mathrm{2}S}\right[{\mathrm{H}}^{+}]/(\mathrm{1}+C_{\mathrm{2}S}\left[{\mathrm{H}}^{+}\right])+[{\mathrm{H}}^{+}]\\ & -\left[{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}\right(\mathrm{2}{C}_{\mathrm{1}}/\left[{\mathrm{H}}^{+}\right]+\mathrm{1})/(\mathrm{1}+C_{\mathrm{1}}/\left[{\mathrm{H}}^{+}\right])=\mathrm{0},\end{array}\end{array}$$

where $C_{\mathrm{1}}=K_{\mathrm{a},{\mathrm{HSO}}_{\mathrm{4}}}\frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{{\mathit{\gamma }}_{\mathrm{2}\mathrm{H}\text{-}{\mathrm{SO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{3}}}$, while $C_{\mathrm{2}S}=\frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{K_{\mathrm{ag},{\mathrm{NH}}_{\mathrm{3}}}{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}$. The only unknown is thus [H⁺], which can thus be solved by bisection solution processes. As shown in the equation, activity coefficients that matter in solving this system include $\frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{{\mathit{\gamma }}_{\mathrm{2}\mathrm{H}\text{-}{\mathrm{SO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{3}}}$ in C₁ and $\frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{HSO}}_{\mathrm{4}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}$ in C_2S.

Calculation principles for subcase D3 (an NH₃–H₂SO₄–HNO₃–H₂O aerosol system). For the subcase D3, the major equilibriums considered is the gas–particle partitioning of ammonia and nitrates of

$$\begin{array}{}\text{(D3)}& {\displaystyle }\begin{array}{rl}& {\mathrm{NH}}_{\mathrm{4}}^{+}\rightleftharpoons {\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{aq}\right)+{\mathrm{H}}^{+},\\ & K_{\mathrm{ag},{\mathrm{NH}}_{\mathrm{3}}}=K_{\mathrm{w}}/K_{\mathrm{bg},{\mathrm{NH}}_{\mathrm{3}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}\left(\mathrm{aq}\right)\left]\right[{\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]{\mathit{\gamma }}_{{\mathrm{H}}^{+}\left(\mathrm{aq}\right)}}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}}{\displaystyle \frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{NO}}_{\mathrm{3}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}},\end{array}\text{(D4)}& {\displaystyle }\begin{array}{rl}& {\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)\rightleftharpoons {\mathrm{H}}^{+}+{\mathrm{NO}}_{\mathrm{3}}^{-},\\ & K_{\mathrm{ag},{\mathrm{HNO}}_{\mathrm{3}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)]{\mathit{\gamma }}_{{\mathrm{H}}^{+}\left(\mathrm{aq}\right)}{\mathit{\gamma }}_{{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)}}{{[\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)]{\mathit{\gamma }}_{{\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)}}}.\end{array}\end{array}$$

Note that in subcase D3, ${\mathit{\gamma }}_{{\mathrm{H}}^{+}}/{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}^{+}}$ is estimated by (${\mathit{\gamma }}_{\mathrm{H}-{\mathrm{NO}}_{\mathrm{3}}}/{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}}{)}^{\mathrm{2}}$, not (${\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{HSO}}_{\mathrm{4}}}/{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}-{\mathrm{HSO}}_{\mathrm{4}}}{)}^{\mathrm{2}}$ as in subcase A2.

These two equilibriums are further combined to be

$$\begin{array}{}\text{(D5a)}& C_{\mathrm{3}}={\displaystyle \frac{{[\mathrm{NO}}_{\mathrm{3}}^{-}{\left(\mathrm{aq}\right)\left]\right[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}{{[\mathrm{HNO}}_{\mathrm{3}}{\left(\mathrm{g}\right)\left]\right[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}}={\displaystyle \frac{K_{\mathrm{ag},{\mathrm{HNO}}_{\mathrm{3}}}}{K_{\mathrm{ag},{\mathrm{NH}}_{\mathrm{3}}}{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}}^{\mathrm{2}}}}.\end{array}$$

As to the charge balance, here only major species are regarded as

$$\begin{array}{}\text{(D5b)}& \begin{array}{rl}\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]& =\mathrm{2}\left[{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}\left(\mathrm{aq}\right)\right]+\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)\right]\\ & =\mathrm{2}[{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}+[{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)].\end{array}\end{array}$$

Combining Eqs. (D5a) and (D5b), at given total nitrate ([HNO₃]_t, namely [NO${}_{\mathrm{3}}^{-}$(aq)] + [HNO₃(g)]) and [NH₃]_t (= [NH${}_{\mathrm{4}}^{+}$(aq)] + [NH₃(g)]) levels, the solution function can be expressed as

$$\begin{array}{}\text{(D6)}& \begin{array}{rl}& {\displaystyle \frac{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}{{[\mathrm{NH}}_{\mathrm{3}}\left(\mathrm{g}\right)]}}{\displaystyle \frac{{[\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)]}{{[\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)]}}-C_{\mathrm{3}}=\mathrm{0},\\ & {\displaystyle \frac{{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]}{{\left(\right[\mathrm{NH}}_{\mathrm{3}}{]}_{\mathrm{t}}-{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\left]\right)}}\\ & {\displaystyle \frac{{\left(\right[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)]-\mathrm{2}{[\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}})}{{\left(\right[\mathrm{HNO}}_{\mathrm{3}}{]}_{\mathrm{t}}+\mathrm{2}{[\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}-{[\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\left]\right)}}-C_{\mathrm{3}}=\mathrm{0},\end{array}\end{array}$$

where the only unknown is [NH${}_{\mathrm{4}}^{+}$(aq)] and can be solved through the bisection method. As shown in the equation, the only activity coefficient that matters in solving this system is (${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}}{)}^{\mathrm{2}}$ in C₃.

Calculation principles for subcase G5 (a Na⁺–NH₃–H₂SO₄–HNO₃–HCl–H₂O aerosol system). For the subcase G5, the major equilibriums considered are the gas–particle partitioning of NH₃, HNO₃ and HCl, while sulfate is considered to exist mainly as [SO${}_{\mathrm{4}}^{\mathrm{2}-}$(aq)]. General derivation processes are similar with D3 and are also detailed in a previous study (Song et al., 2018). Briefly, the key equilibriums include that of HNO₃ (Eq. D4) and HCl of

$$\begin{array}{}\text{(D7)}& \begin{array}{rl}& \mathrm{HCl}\left(\mathrm{g}\right)\rightleftharpoons {\mathrm{H}}^{+}\left(\mathrm{aq}\right)+{\mathrm{Cl}}^{-}\left(\mathrm{aq}\right),\\ & K_{\mathrm{ag},\mathrm{HCl}}={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{Cl}}^{-}\left(\mathrm{aq}\right)]{\mathit{\gamma }}_{{\mathrm{H}}^{+}\left(\mathrm{aq}\right)}{\mathit{\gamma }}_{{\mathrm{Cl}}^{-}\left(\mathrm{aq}\right)}}{\left[\mathrm{HCl}\right(\mathrm{g}\left)\right]{\mathit{\gamma }}_{\mathrm{HCl}\left(\mathrm{g}\right)}}}\\ & ={\displaystyle \frac{{[\mathrm{H}}^{+}{\left(\mathrm{aq}\right)\left]\right[\mathrm{Cl}}^{-}\left(\mathrm{aq}\right)]{\mathit{\gamma }}_{\mathrm{H}\text{-}\mathrm{Cl}}^{\mathrm{2}}}{\left[\mathrm{HCl}\right(\mathrm{g}\left)\right]}},\end{array}\end{array}$$

which can be combined into

$$\begin{array}{}\text{(D8)}& \begin{array}{rl}{C}_{\mathrm{4}}& ={\mathit{\xi }}_{{\mathrm{HNO}}_{\mathrm{3}}}/{\mathit{\xi }}_{\mathrm{HCl}}={\displaystyle \frac{{[\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)\left]\right[\mathrm{HCl}\left(\mathrm{g}\right)]}{{[\mathrm{Cl}}^{-}{\left(\mathrm{aq}\right)\left]\right[\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)]}}\\ & ={\displaystyle \frac{K_{\mathrm{ag},{\mathrm{HNO}}_{\mathrm{3}}}}{K_{\mathrm{ag},\mathrm{HCl}}}}{\displaystyle \frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}\mathrm{Cl}}^{\mathrm{2}}}{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{NO}}_{\mathrm{3}}}^{\mathrm{2}}}}.\end{array}\end{array}$$

Therefore [NO${}_{\mathrm{3}}^{-}$ (aq)] and [HNO₃(g)] (=[HNO₃]_t − [NO${}_{\mathrm{3}}^{-}$ (aq)]) can be solved at known assumed [Cl⁻(aq)].

The [NH${}_{\mathrm{4}}^{+}$(aq)] associated with Cl⁻(aq) and NO${}_{\mathrm{3}}^{-}$(aq), [NH${}_{\mathrm{4}}^{+}$(aq)]_NC, is solved by

$$\begin{array}{}\text{(D9)}& \left(\right[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right){]}_{\mathrm{NC}}{)}^{\mathrm{2}}-\mathrm{B}[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right){]}_{\mathrm{NC}}+C=\mathrm{0},\end{array}$$

where

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }B=[{\mathrm{NH}}_{\mathrm{3}}{]}_{\mathrm{t}}+[{\mathrm{Na}}^{+}]-\mathrm{2}[{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}+\left[{\mathrm{Cl}}^{-}\left(\mathrm{aq}\right)\right]\\ {\displaystyle }& {\displaystyle }+\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)\right]+C_{\mathrm{2}\mathrm{N}}^{-\mathrm{1}},\\ {\displaystyle }& {\displaystyle }C=\left(\right[{\mathrm{NH}}_{\mathrm{3}}{]}_{\mathrm{t}}+\left[{\mathrm{Na}}^{+}\right]-\mathrm{2}\left[{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}\right)\left(\right[{\mathrm{Cl}}^{-}\left(\mathrm{aq}\right)]\\ {\displaystyle }& {\displaystyle }+\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\left(\mathrm{aq}\right)\right])-C_{\mathrm{2}\mathrm{N}}^{-\mathrm{1}}(\mathrm{2}[{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}-[{\mathrm{Na}}^{+}\left]\right),\end{array}$$

where $C_{\mathrm{2}\mathrm{N}}=\frac{{\mathit{\gamma }}_{\mathrm{H}\text{-}{\mathrm{NO}}_{\mathrm{3}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}{K_{\mathrm{ag},{\mathrm{NH}}_{\mathrm{3}}}{\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}\left(\mathrm{aq}\right)}^{\mathrm{2}}}$. And with [NH${}_{\mathrm{4}}^{+}$(aq)]_NC, we have

$$\begin{array}{c}{\displaystyle }\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)\right]=[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right){]}_{\mathrm{NC}}+\mathrm{2}[{\mathrm{H}}_{\mathrm{2}}{\mathrm{SO}}_{\mathrm{4}}{]}_{\mathrm{t}}-\left[{\mathrm{Na}}^{+}\right],\\ {\displaystyle }\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{g}\right)\right]=[{\mathrm{NH}}_{\mathrm{3}}{]}_{\mathrm{t}}-[{\mathrm{NH}}_{\mathrm{4}}^{+}\left(\mathrm{aq}\right)].\end{array}$$

The system then solves the equation sets through the bisection method by assuming a series of [Cl⁻(aq)] levels.

As shown in the equations above, activity coefficients that matter in solving this system (Eqs. D8–D9) include ${\mathit{\gamma }}_{\mathrm{H}-{\mathrm{NO}}_{\mathrm{3}}}$, γ_H−Cl and ${\mathit{\gamma }}_{{\mathrm{NH}}_{\mathrm{4}}\text{-}{\mathrm{NO}}_{\mathrm{3}}}$.

Appendix E: Information for the assumed calibration year of 2012 in Toronto site (Fig. E1).

Figure E1 Monthly variation in (a, c) NO₃- fraction in anions(aq) and (b, d) the corresponding non-ideality correction factors for the Toronto site in 2012. The data are estimated by (a–b) the E-AIM model and (c–d) the ISORROPIA model. The black dash lines represent the annual mean levels. The box and whiskers represent the 10th, 25th, 50th, 75th and 90th percentiles, while the red markers represent the monthly means.