Molecular-level nucleation mechanism of iodic acid and methanesulfonic acid

Both iodic acid (HIO₃, IA) and methanesulfonic acid (CH₃S(O)₂OH, MSA) have been identified by field studies as important precursors of new particle formation (NPF) in marine areas. However, the mechanism of NPF in which IA and MSA are jointly involved is still unclear. Hence, we investigated the IA-MSA nucleation system under different atmospheric conditions and uncovered the corresponding nucleating mechanism at the molecular level for the first time, using a quantum chemical approach and Atmospheric Cluster Dynamics Code (ACDC). The findings show that the pure-IA nucleation rate was much lower than the results of CLOUD (Cosmics Leaving Outdoor Droplets) experiments. MSA can promote IA cluster formation through stabilizing IA via both hydrogen and halogen bonds, especially under conditions with lower temperatures, sparse IA, and rich MSA. However, the nucleation rate of the IA-MSA mechanism is much lower than that of field observations, indicating that the effect of additional nucleation precursors needs to be considered (e.g., H₂SO₄, HIO₂, NH₃, and amines). The IA-MSA nucleation mechanism revealed in this study may help to gain insight into the joint effect of marine sulfur- and iodine-containing components on marine NPF.

1 Introduction

Marine aerosols, being the primary natural aerosol system (O'Dowd and De Leeuw, 2007), significantly affect global radiation balance and climate by regulating cloud properties as cloud condensation nuclei (CCN) (Takegawa et al., 2020; IPCC, 2013). Nearly half of the CCN originate from new particle formation (NPF) via the gas-to-particle conversion (Merikanto et al., 2009; Yu and Luo, 2009). As a major source of CCN globally, NPF mainly consists of the nucleation of gaseous molecules and the subsequent growth of the formed clusters (Kulmala et al., 2013; Kulmala, 2003; Zhang, 2010). Although extensive studies have provided observational evidence of NPF events in the coastal zone, open ocean, and even ice-covered polar regions (Zheng et al., 2021; Sipilä et al., 2016; Yu et al., 2019; Baccarini et al., 2020), the corresponding NPF mechanisms at the molecular level remain poorly understood, stemming from the lack of chemical speciation in the initial nucleating steps.

Marine NPF, particularly in remote marine areas, is more affected by biological emissions compared to inland ones with anthropogenic influence (Kerminen et al., 2018). Historically, sulfur-containing species originating from ocean-emitted dimethyl sulfide (DMS) have long been identified as significant components of marine aerosols (Charlson et al., 1987; Shaw, 1983; Bates et al., 1992). Methanesulfonic acid (CH₃S(O)₂OH, MSA), as a well-known oxidation product of DMS (Chen et al., 2018; Hatakeyama et al., 1982), is widely dispersed throughout the world's oceans and has considerable atmospheric concentrations (Chen et al., 2018), comparable to or higher than sulfuric acid (SA), i.e., [MSA] $/$ [SA] = 10 %–250 % (Berresheim et al., 2002; Davis et al., 1998; Eisele and Tanner, 1993). Moreover, MSA has been experimentally demonstrated to be a significant nucleating precursor in coastal and remote oceans (Dawson et al., 2012; Hodshire et al., 2019; Karl et al., 2007). Along with stricter global controls on anthropogenic SO₂ emissions, the impact of MSA on NPF will become increasingly significant in the future (Perraud et al., 2015), particularly in marine areas.

In addition to the above sulfur precursors, recent experimental and theoretical studies (He et al., 2021; Martín et al., 2020; Xia et al., 2020; Rong et al., 2020) have also recognized the critical role of iodine compounds in marine NPF processes. According to the field studies (O'Dowd et al., 2002; Sipilä et al., 2016), the observed intense NPF events occur during low tide and are accompanied by a significant increase in iodic acid (HIO₃, IA) concentration in the coastal Mace Head, Ireland, indicating that the coastal NPF is primarily driven by subsequential addition of IA and involves the participation of I₂O₅. More recently, He et al. (2021) demonstrated experimentally that, in addition to IA and I₂O₅, iodous acid (HIO₂) and I₂O₄ are also involved in the cluster formation process, with HIO₂ playing a key role in the stabilization of neutral IA clusters. Also, recent evidence suggests that the NPF events in the sea-ice-covered Arctic region are also mainly driven by IA (Baccarini et al., 2020). Notably, in addition to IA, significant concentrations of MSA were observed during marine NPF events (Beck et al., 2021). Although MSA and IA were detected in the smallest clusters (Beck et al., 2021), it is still unknown whether they could be simultaneously involved in the early nucleation process. If so, their joint nucleation mechanism and the corresponding regions affected by that mechanism need to be further elucidated.

Herein, the high-level quantum chemical calculations combined with Atmospheric Cluster Dynamics Code (ACDC) (McGrath et al., 2012) were employed to simulate the nucleating process of the (IA)_x ⋅ (MSA)_y system (where $\mathrm{0}\le x\le \mathrm{6}$, $\mathrm{0}\le y\le \mathrm{3}$, $\mathrm{1}<x+y\le \mathrm{6}$). Under different atmospheric conditions (temperature and precursor concentration), a series of ACDC simulations were carried out to explore (i) the binding nature of IA and MSA, (ii) the joint effects of IA and MSA on the nucleation, and (iii) which conditions are more affected by the IA-MSA mechanism. The current work may contribute to developing a more comprehensive marine NPF mechanism and explaining some missing sources of particles in marine environments.

2 Methods

2.1 Quantum chemistry calculation

All structure optimizations with tight convergence criteria and frequency calculations with density functional theory (DFT) were carried out by the Gaussian 09 package (Frisch et al., 2009). Considering the variety of possible isomers of multimolecular clusters, a systematic multistep conformation search was employed here to locate the lowest-energy cluster structures. The structures of the pure-IA clusters employed in this study refer to those in the study by Rong et al. (2020). A new (IA)₆ cluster with lower energy was found here. For each studied IA-MSA cluster, the artificial bee algorithm combining the UFF force field (Rappé et al., 1992) was adopted to yield 1000 initial configurations from 5000 generations by ABCluster software (Zhang and Dolg, 2015). After pre-optimization by the PM7 semiempirical method (Stewart, 2013) with MOPAC2016 (Stewart, 2016), the 100 structures with lower energies were left for further optimization at the ωB97X-D/6-31 + G^* + Lanl2DZ (for iodine) level of theory due to the best performance of the ωB97X-D functional in studying atmospheric clusters (Elm and Kristensen, 2017; Schmitz and Elm, 2020). The final global minima were reoptimized by ωB97X-D functional with the 6-311 $++$ G(3df,3pd) (Francl et al., 1982) basis set for H, C, O, and S atoms and aug-cc-pVTZ-PP with ECP28MDF for I atom (Peterson et al., 2003), and they were identified with the lowest Gibbs formation free energy (ΔG). It is noteworthy that the larger aug-cc-pVTZ-PP basis set (for iodine) was employed in the present study compared to the aug-cc-pVDZ-PP basis set (for iodine) in our previous work (Rong et al., 2020), because higher level of theory usually implies a better calculation accuracy.

The single-point correction was further performed by the RI-CC2 method (Hattig and Weigend, 2000) with aug-cc-pVTZ (for H, C, O) + aug-cc-pV(T+d)Z (for S) + aug-cc-pVTZ-PP with ECP28MDF (for I) basis set using the TURBOMOLE program (Dunning et al., 2001; Ahlrichs et al., 1989), since the ACDC simulations based on RI-CC2 values are in good agreement with the experimental results (Lu et al., 2020; Kürten et al., 2018; Li et al., 2020; Almeida et al., 2013). Herein, the Gibbs formation free energies (ΔG, kcal mol⁻¹) of the studied clusters were calculated as Eq. (1):

$$\begin{array}{}\text{(1)}& \mathrm{\Delta }G=\mathrm{\Delta }{E}_{\text{RI-CC2}}+\mathrm{\Delta }{G}_{\mathrm{thermal}}^{\mathit{\omega }\text{B97X-D}},\end{array}$$

where ΔE_RI-CC2 is the electronic contribution and $\mathrm{\Delta }{G}_{\mathrm{thermal}}^{\mathit{\omega }\text{B97X-D}}$ is the thermal contribution to Gibbs free energy. For subsequent clustering kinetic simulations at different temperatures, the ΔG values of clusters ranging from 218 to 298 K were calculated by Shermo 2.0 (Lu and Chen, 2021) and collected in Table S3.

2.2 Wave function analysis

To better understand the interactions between IA and MSA, the bonding nature was investigated through wave function analysis using Multiwfn 3.7 (Lu and Chen, 2012). Specifically, the molecular electrostatic potential (ESP) was calculated for IA and MSA, which facilitates understanding their potential interaction sites. Moreover, natural bond orbital (NBO) analysis (Reed et al., 1988) was carried out to give a detailed insight into intermolecular interactions. Based on the final identified stable clusters, the NBO information calculated by Gaussian 09 is resolved by Multiwfn, and the key interactive orbitals are visualized by VMD 1.9.3 (Humphrey et al., 1996). To further quantify the binding strength, electron density ρ(r), Laplacian electron density ∇²ρ(r), energy density H(r) at corresponding bond critical points (BCPs) based on atoms in molecules (AIM) theory (Becke, 2007; Lane et al., 2013) were also calculated in this work (Table S1 in the Supplement).

2.3 Atmospheric cluster dynamic simulations

Simulation of the nucleation process of the IA-MSA system is achieved by the Atmospheric Clusters Dynamic Code (ACDC) (McGrath et al., 2012). Specifically, the ACDC derives the steady-state concentration and cluster formation rates by solving the birth–death equations (Eq. 2).

$$\begin{array}{}\text{(2)}& \begin{array}{rl}{\displaystyle \frac{\mathrm{d}{c}_i}{\mathrm{d}t}}& ={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\sum _{j<i}{\mathit{\beta }}_{j,\left(i-j\right)}{C}_j{C}_{\left(i-j\right)}+\sum _j{\mathit{\gamma }}_{\left(i+j\right)\to i}{C}_{i+j}\\ & -\sum _j{\mathit{\beta }}_{i,j}{C}_i{C}_j-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\sum _{j<i}{\mathit{\gamma }}_{i\to j}{C}_i+Q_i-S_i,\end{array}\end{array}$$

where C_i refers to the cluster i concentration, β_i,j is the collision rate coefficient between clusters i and j, γ_i→j is the evaporation rate coefficient of smaller cluster j from the parent cluster i, and Q_i and S_i are the outside source and loss term of cluster i, respectively. β_i,j is calculated based on the kinetic gas theory, which is given as

$$\begin{array}{}\text{(3)}& {\mathit{\beta }}_{i,j}={\left({\displaystyle \frac{\mathrm{3}}{\mathrm{4}\mathit{\pi }}}\right)}^{\mathrm{1}/\mathrm{6}}{\left({\displaystyle \frac{\mathrm{6}{k}_{\mathrm{B}}T}{m_i}}+{\displaystyle \frac{\mathrm{6}{k}_{\mathrm{B}}T}{m_j}}\right)}^{\mathrm{1}/\mathrm{2}}{\left(V_i^{\mathrm{1}/\mathrm{3}}+V_j^{\mathrm{1}/\mathrm{3}}\right)}^{\mathrm{2}},\end{array}$$

where V_i and m_i are the volume and mass of cluster i, respectively. k_B is the Boltzmann constant, and T is the temperature. Equation (3) is derived from the hard-sphere collision theory where $V_i=\mathrm{3}/\mathrm{4}\times \mathit{\pi }\times (D_i/\mathrm{2}{)}^{\mathrm{3}}$. The diameter D_i of cluster i is calculated by Multiwfn (Lu and Chen, 2012). Evaporation rate coefficients, ${\mathit{\gamma }}_{(i+j)\to i,j}$, are derived from ΔG of clusters and the corresponding collision rate coefficients based on the detailed balance assumption (McGrath et al., 2012):

$$\begin{array}{}\text{(4)}& {\mathit{\gamma }}_{\left(i+j\right)\to i,j}={\mathit{\beta }}_{i,j}{\displaystyle \frac{P_{\mathrm{ref}}}{k_{\mathrm{B}}T}}\exp \left({\displaystyle \frac{\mathrm{\Delta }{G}_{i+j}-\mathrm{\Delta }{G}_i-\mathrm{\Delta }{G}_j}{k_{\mathrm{B}}T}}\right),\end{array}$$

where P_ref is the reference pressure (1 atm) at which the Gibbs free energies were determined, and ΔG_i is the Gibbs formation free energy of the formation of cluster i from the corresponding monomers.

In the present study, the ACDC simulations only modeled the neutral cluster formation process and did not consider the charge nor the effect of water. Since IA is weakly bound to water, it is less inclined to exist as hydration of IA in tropospheric conditions (Khanniche et al., 2016). Meanwhile, the nucleation efficiency of MSA and water is low (Arquero et al., 2017). Thus, the effect of water on the conclusion in this study is limited. The settings of the boundary conditions of the ACDC simulations are discussed in Sect. S1 (Supplement) and summarized in Table S5.

3 Results

3.1 Cluster conformational analysis

The obtained most stable structures of (IA)_x ⋅ (MSA)_y ($\mathrm{0}\le x\le \mathrm{6}$, $\mathrm{0}\le y\le \mathrm{3}$, $\mathrm{1}<x+y\le \mathrm{6}$) clusters are presented in Fig. S1 in the Supplement, and the corresponding Cartesian coordinates are collected in Table S7 in the Supplement. To investigate the intermolecular bonding potential of IA and MSA, the electrostatic potential (ESP) was calculated to analyze their potential interaction sites.

As shown in Fig. 1a, IA has positive ESPs (red region) surrounding its −OH group, with a maximum value of +59.04 kcal mol⁻¹, making the −OH group an effective hydrogen bond (HB) donor. And IA's two terminal oxygens with negative ESPs (−29.09 and −29.47 kcal mol⁻¹) can serve as HB acceptors. Similarly, the −OH group of MSA has the strongest electrophilicity (ESP value of +63.86 kcal mol⁻¹) as the HB donor, while its terminal O atom has strong nucleophilicity as the HB acceptor, due to its lone pair of electrons. In this case, IA and MSA can directly bind with each other via HBs. Moreover, IA possesses positive charge localization (the so-called δ hole) with a maximal ESP value of +51.87 kcal mol⁻¹ at the end of the iodine atom along the O–I direction. This electron-deficient region tends to attract the electron-rich oxygen atom of the MSA to form the halogen bonds (XBs) O–I^…O (green band line in Fig. 1a). From the skeletal formula presented in Fig. 1a, we can know that a formed (IA)₁ ⋅ (MSA)₁ cluster is stabilized by both HBs and XBs. A similar situation has also been found in the larger IA-MSA clusters in Fig. S1.

Figure 1 (a) The ESP-mapped molecular vdW (van der Waals) surfaces of iodic acid (IA) and methanesulfonic acid (MSA). The pink, red, yellow, cyan, and white spheres represent I, O, S, C, and H atoms, respectively. The yellow and cyan dots indicate the positions of maximums and minimums of ESP (in kcal mol⁻¹), respectively. (b) The donor–acceptor NBOs involved in the (IA)₁ ⋅ (MSA)₁ cluster. LP indicates the lone-pair orbitals, and δ^∗ indicates the antibonding orbitals.

Additionally, as illustrated in Fig. 1b, natural bond orbital (NBO) analysis was performed to reveal the bonding nature of IA and MSA. For the formed O–I^…O halogen bond, the lone-pair orbital LP(O) in the terminal oxygen atom of MSA acts as an electron donor, while the antibonding orbital δ^∗(O–I) in IA is the electron acceptor. Essentially, halogen bonding originates from the interactions between LP(O) and δ^∗(O–I) orbitals, accompanied by intermolecular charge transferring from LP(O) to δ^∗(O–I). In the case of the O–H^…O hydrogen bond, the LP(O) of IA serves as the donor orbital, and δ^∗(O–H) of MSA is the acceptor orbital; the charge shifts from LP(O) to δ^∗(O–H). The ESP and NBO results indicate that IA and MSA are capable of forming both HBs and XBs and have the potential to form stable clusters.

To quantify the bonding strength of HBs or XBs within the studied IA-MSA clusters (Fig. S1), the bonding properties, including electron density ρ(r), Laplacian electron density ∇²ρ(r), and energy density H(r) at the bond critical points (BCPs), are calculated based on the AIM methodology (Becke, 2007; Lane et al., 2013) and collected in Table S1. For O–I^…O XBs, the ρ(r), ∇²ρ(r), and H(r) values at the BCPs are in the ranges of 0.0143 to 0.0849, 0.0409 to 0.1589, and −0.0265 to 0.0019 a.u., respectively. As for O–H^…O HBs, the ρ(r), ∇²ρ(r), and H(r) values at the BCPs are in the ranges of 0.0178 to 0.0796, 0.0615 to 0.1141, and −0.0032 to −0.0332 a.u., respectively. The electron density ρ(r) is generally positively correlated with the bond strength, and the ρ(r) values (Table S1) are well within the specified ρ(r) range of HBs (0.002 to 0.040 a.u.) (Grabowski, 2004; Koch and Popelier, 1995), indicating that all the O–H^…O non-covalent interactions are indeed HBs. Moreover, according to the classification of HBs (Rozas et al., 2000), all the HBs formed within IA-MSA clusters are medium HBs (12.0<E (interaction energy) <24.0 kcal mol⁻¹) with ∇²ρ(r)>0 and H(r)<0. Overall, the conformational analysis suggests that MSA can stabilize IA clusters by forming relatively strong non-covalent interactions such as HBs and XBs and thus has the potential to form relatively stable clusters with IA.

3.2 Cluster stability analysis

To evaluate the thermodynamic stability of formed IA-MSA clusters, the Gibbs formation free energy (ΔG, kcal mol⁻¹) of each studied (IA)_x ⋅ (MSA)_y ($\mathrm{0}\le x\le \mathrm{6}$, $\mathrm{0}\le y\le \mathrm{3}$, $\mathrm{1}<x+y\le \mathrm{6}$) cluster at T = 218–298 K ranging from boundary layer to free troposphere (Williamson et al., 2019) and p=1 atm was calculated by Eq. (1) and shown in Table S3. In the present study, the analysis and discussion of the simulation results are mainly at T=278 K.

As shown in Fig. 2a, the ΔG values of IA-MSA clusters at 278 K decrease with increasing cluster size, indicating that the cluster growth process is energetically favorable. And the same trend is also observed at 298 K (Fig. S2) and 258 K (Fig. S3). The ΔG of the (IA)_x ⋅ (MSA)₁ (x = 1–5) clusters, are 7.71–15.67 kcal mol⁻¹ lower than those of the corresponding (IA)_x clusters, indicating that pure-IA clusters could potentially grow by binding with MSA. Moreover, the corresponding total evaporation rate coefficients (∑γ, s⁻¹) of clusters were calculated at 278 K by Eq. (4) and presented in Fig. 2b and Table S4. In general, a lower ∑γ value indicates greater cluster stability. As shown in Fig. 2b, the ∑γ values of larger clusters, (IA)_4–6 and (IA)_3–4 ⋅ (MSA)₂, are significantly lower than those of the corresponding initial, small-sized clusters, indicating that the stability increases during the cluster growth. Considering the competition between collision and evaporation during the clustering process, the ratio of collision frequencies versus total evaporation rate coefficients ($\mathit{\beta }C/\sum \mathit{\gamma }$) was calculated to access the probability of cluster growth. If $\mathit{\beta }C/\sum \mathit{\gamma }>\mathrm{1}$, further growth of cluster by colliding IA or MSA molecules can be assumed to dominate over cluster evaporation (details in the Sect. S1 of the Supplement). Figure 2c presents the minimum values of ${\mathit{\beta }}_{\mathrm{IA}}{C}_{\mathrm{IA}}/\sum \mathit{\gamma }$ for the studied clusters at lowest C_IA (10⁶ molec. cm⁻³), where β_IA is the rate coefficient of cluster collision with IA monomer, and C_IA is the concentration of IA monomer. Figure 2d presents the results of the collision with MSA (${\mathit{\beta }}_{\mathrm{MSA}}{C}_{\mathrm{MSA}}/\sum \mathit{\gamma }$). Similarly, the minimum values of ${\mathit{\beta }}_{\mathrm{IA}}{C}_{\mathrm{IA}}/\sum \mathit{\gamma }$ at 298 and 258 K are presented in Figs. S2 and S3, respectively. Among these clusters, the largest (IA)₄ ⋅ (MSA)₂ and (IA)₆ clusters ($\mathit{\beta }C/\sum \mathit{\gamma }>\mathrm{1}$) incline to collide with IA monomer (or MSA monomer) to grow out of the simulated system. As a result, the fluxes for clusters with larger size than (IA)₄ ⋅ (MSA)₂ and (IA)₆ are counted in the cluster formation rate J.

Figure 2 (a) Gibbs formation free energy (ΔG, kcal mol⁻¹) of the (IA)_x ⋅ (MSA)_y ($\mathrm{0}\le x\le \mathrm{6}$, $\mathrm{0}\le y\le \mathrm{3}$, $\mathrm{1}<x+y\le \mathrm{6}$) clusters calculated at the RI-CC2/aug-cc-pV(T+d)Z(-PP)//ωB97X-D/6-311 $++$ G(3df,3pd) + aug-cc-pVTZ-PP level of theory; T=278 K, and p=1 atm. (b) The total evaporation rate coefficients (∑γ, s⁻¹) and the ratios of collision frequencies with (c) IA monomer or (d) MSA monomer versus total evaporation rate coefficients (${\mathit{\beta }}_{\mathrm{IA}}{C}_{\mathrm{IA}}/\sum \mathit{\gamma }$ or ${\mathit{\beta }}_{\mathrm{MSA}}{C}_{\mathrm{MSA}}/\sum \mathit{\gamma }$) of the corresponding clusters. β is the rate coefficient of cluster collision with monomer, and C is the monomer concentration.

3.3 Cluster formation rates

To comprehensively explore the effect of MSA on IA cluster formation kinetically, the IA-MSA cluster formation rate J (cm⁻³ s⁻¹) was simulated under different atmospheric conditions using ACDC. Herein, we first explored the changes of J after the intervention of different concentrations of MSA ([MSA]), using the pure-IA system as a reference. Based on the field measurement, [IA] in the ACDC simulation is set to be in the range of 10⁶–10⁸ molec. cm⁻³ (Sipilä et al., 2016); [MSA] = 10⁶, 10⁷, and 10⁸ molec. cm⁻³ (Chen et al., 2018; Berresheim et al., 2002; Davis et al., 1998). As for the setting of condensation sink coefficient (CS), the different CS values ($\mathrm{1.0}\times {\mathrm{10}}^{-\mathrm{4}}$–2.6 × 10⁻³ s⁻¹) have an impact on the simulated J, especially in the case of low J (Fig. S4), but less on presenting the promotion of MSA on IA cluster formation and the main conclusions of this study. Hence, the CS is chosen as a typical coastal value ($\mathrm{2.0}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹) (Dal Maso et al., 2002), which is uniform for all clusters.

Figure 3 Simulated cluster formation rates J (cm⁻³ s⁻¹) as a function of iodic acid concentration [IA], with different concentrations of methanesulfonic acid [MSA] of 10⁶ (blue), 10⁷ (yellow), 10⁸ (red), and 0 molec. cm⁻³ (purple, “pure-IA”), at T = 278 K, and CS = $\mathrm{2.0}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹. The gray diamonds are the measured rate data from the CLOUD experiment at T=283 K (He et al., 2021).

As shown in Fig. 3, the simulated J of pure-IA nucleation (purple line) is much lower than the rate obtained from the CLOUD experiment (He et al., 2021). The J values of the IA-MSA system with varying [MSA] (red, yellow, and blue lines) are all higher than those of the pure-IA system (purple line). Particularly, at a lower [IA] of 10⁶ molec. cm⁻³, the involvement of MSA results in a greater boost in J. Briefly, MSA can promote J of IA clusters to a higher level, which is a reflection of the stabilizing effect of MSA on IA clusters. However, the J of IA-MSA nucleation was still much less than the experimental results (He et al., 2021), even at a high [MSA] (10⁸ molec. cm⁻³). The large rate difference suggests that MSA stabilizes IA less efficiently than the potential iodine-containing components.

To quantify such enhancement of MSA on J, here we defined an enhancement strength R as the following Eq. (5):

$$\begin{array}{}\text{(5)}& R={\displaystyle \frac{J_{\text{IA-MSA}}}{J_{\text{pure-IA}}}}={\displaystyle \frac{J\left(\right[\mathrm{IA}]=x,[\mathrm{MSA}]=y)}{J\left(\right[\mathrm{IA}]=x)}},\end{array}$$

where J_IA-MSA and J_pure-IA indicate the J of IA-MSA and pure-IA nucleating systems, respectively. x and y are the atmospheric concentrations of IA and MSA, respectively.

During nucleating processes, variations in ambient conditions (precursor concentration and temperature) can affect J_IA-MSA and J_pure-IA as well as the R of MSA. Herein, the simulations were performed in a wide range of atmospheric temperatures (T = 218 to 298 K) and concentrations of IA (10⁶–10⁸ molec. cm⁻³) and MSA (10⁶ to 10⁸ molec. cm⁻³).

Figure 4 Enhancement strength R of MSA on cluster formation rates under different atmospheric conditions: T = 218 to 298 K (blue cones, [IA] = 10⁷ and [MSA] = 10⁷ molec. cm⁻³), [MSA] = 10⁶ to 10⁸ molec. cm⁻³ (red cones, T = 278 K and [IA] = 10⁷ molec. cm⁻³), [IA] = 10⁶ to 10⁸ molec. cm⁻³ (purple cones, T = 278 K and [MSA] = 10⁷ molec. cm⁻³), and CS = 2.0 × 10⁻³ s⁻¹.

As shown in Fig. 4, the enhancement strength, R, of MSA decreases with increasing [IA] (10⁶ → 10⁸ molec. cm⁻³), under the condition of T = 278 K and [MSA] = 10⁷ molec. cm⁻³ (purple line). The specific R values are summarized in Table S6. This is because the contribution of pure-IA clusters to nucleation becomes higher with increasing [IA], thereby diminishing that of IA-MSA clusters (smaller R). Conversely at lower [IA] (10⁶ molec. cm⁻³), the effect of R of MSA on J could reach 211-fold, even when [MSA] is only at a median value (10⁷ molec. cm⁻³). The R increases with increasing [MSA] (10⁶ → 10⁸ molec. cm⁻³) (orange line) due to more IA-MSA clusters formed. Interestingly, as the temperature decreases from 298 to 218 K (blue line), R first increases (298 → 238 K) and then decreases (238 → 218 K). During the temperature range from 298 to 238 K, the decrease in temperature diminishes cluster evaporation (Eq. 4), which in turn promotes IA-MSA cluster formation and leads to an increase in R. When the temperature is very low, between 218 and 238 K, the effect of cluster evaporation is almost negligible, and the nucleation process is primarily limited by collisions between clusters or molecules, namely, the kinetic limit process. In this case, the lower T reduces the collision rate and thus results in a decrease in R. The numerator, J_IA-MSA, in Eq. (5) is affected by the formation of both pure-IA clusters and IA-MSA clusters, while the denominator, J_pure-IA, is only affected by the formation of pure-IA clusters. When the overall intermolecular collision rate between IA and IA, as well as IA and MSA, is reduced due to the decrease in temperature, the numerator would be affected more than the denominator in Eq. (5), which in turn leads to a reduced R. As a result of the above analysis, the effect of the IA-MSA system on the nucleation process varies with the [IA], [MSA], and T, and this is particularly important in regions with lower T, sparse IA, and rich MSA.

Figure 5 (a) Main cluster growth pathway of the IA-MSA nucleating system at T=278 K, CS = $\mathrm{2.0}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹, [IA] = 10⁷, and [MSA] = 5×10⁶ molec. cm⁻³. The black and orange arrows refer to the pathways of colliding with IA and MSA, respectively, where the dashed arrows indicate the evaporation of MSA. (b) Branch ratio of IA-MSA (orange pieces) and pure-IA (purple pieces) growth pathway under varying [MSA] (10⁶–10⁷ molec. cm⁻³) and [IA] (10⁶–10⁸ molec. cm⁻³).

3.4 Cluster growth pathways

According to the analysis above, MSA can stabilize IA clusters, thereby enhancing cluster formation rate. However, the mechanism of how MSA and IA jointly contribute to cluster formation is still unclear. Thus, the detailed cluster growth pathways were tracked by ACDC and shown in Fig. 5a.

The main clustering pathways can be divided into two types: (i) IA self-nucleation and (ii) IA-MSA cluster formation. The studied clusters that did not appear in the cluster growth pathway are mainly due to their low stability. For the IA self-nucleation pathway, cluster growth proceeds mainly via the collisional binding of IA monomers ((IA)${}_{\mathrm{1}\to \mathrm{2}\to \mathrm{3}\to \mathrm{4}\to \mathrm{5}\to \mathrm{6}}$), which is consistent with the reported pure-IA nucleation mechanism (Rong et al., 2020). For the IA-MSA pathway, it starts from the heterodimer (IA)₁ ⋅ (MSA)₁ and then grows primarily through IA addition, resulting in the (IA)₄ ⋅ (MSA)₂ clusters with sufficient stability to grow out of the simulated system (Fig. 2c). The results suggest that MSA can directly participate in the IA-involved nucleation by forming relatively stable IA-MSA clusters.

Figure 6 The simulated cluster formation rate J (cm⁻³ s⁻¹) of the IA-MSA system at different temperatures (a) 218, (b) 238, (c) 258, and (d) 298 K; [IA] = 10⁶–10⁸ molec. cm⁻³; [MSA] = 0, 10⁶, 10⁷, and 10⁸ molec. cm⁻³; and CS = $\mathrm{2.0}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹.

In the atmosphere, the distribution of IA and MSA varies by region, affecting the contribution of IA-MSA clustering pathways accordingly. Hence, the branch ratios of flux out through the IA-MSA path (orange pieces) and pure-IA path (purple pieces) at varying [MSA] (10⁶–10⁷ molec. cm⁻³) and [IA] (10⁶–10⁸ molec. cm⁻³) are presented in Fig. 5b to access the impact of the IA-MSA mechanism. As shown in Fig. 5b, the branch ratio of IA-MSA and pure-IA path is highly dependent on [MSA] and [IA]. At the conditions of T = 278 K, CS = $\mathrm{2.0}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹, and [IA] = 10⁷ molec. cm⁻³, the contribution of the IA-MSA path increases from 1 % to 66 % with increasing [MSA]. Additionally, given the uneven distribution of IA, the analysis was further carried out within the atmospherically relevant range of [IA] (10⁶–10⁸ molec. cm⁻³). The results show that the contribution of the IA-MSA path decreases from 94 % to 2 % with increasing [IA] (10⁶–10⁸ molec. cm⁻³). These findings indicate that the IA-MSA mechanism contributes more in regions with higher [MSA] and lower [IA]. Furthermore, the branch ratio was calculated based on field conditions (temperatures and [IA]) reported by He et al. (2021) and presented in Fig. S5. The results indicate that the IA-MSA mechanism does have stronger effects in polar regions than in mid-latitude coastal regions due to lower temperatures, which is also consistent with the above findings.

Most of the analysis above was performed at 278 K. To further probe the impact of temperature on J systematically, Fig. 6 presents the simulated J at additional temperatures (218, 238, 258 and 298 K), [IA] = 10⁶–10⁸ molec. cm⁻³, [MSA] = 10⁶ (red line), 10⁷ (yellow line), and 10⁸ molec. cm⁻³ (purple line). At a relatively high T=298 K (Fig. 6d), the improvement on J by the addition of MSA was not significant compared to the pure-IA system, except at higher [MSA] = 10⁸ molec. cm⁻³ and relatively lower [IA]. At lower T=258 K (Fig. 6c), the enhancement in J by MSA is stronger in all cases except at lowest [MSA] = 10⁶ molec. cm⁻³. Moreover, such a boost in J was further enhanced at 238 K (Fig. 6b). Lower concentrations of MSA (10⁶ molec. cm⁻³) also significantly promote the formation of IA clusters, mainly because the low temperature weakens the cluster evaporation.

Interestingly, the comparison of the simulations at 218 K (Fig. 6a) and 238 K (Fig. 6b) shows that the decrease in temperature does not further improve J to a higher level. Instead, there is a decrease of the enhancement of MSA at higher [IA] = 10⁷–10⁸ molec. cm⁻³. The reason for this phenomenon, which was also discussed in Sect. 3.3, is that at such low temperatures cluster growth is more dependent on collisions and less on evaporation. And the low temperature reduces the rate of collision between clusters or molecules.

Compared to the field observations at Mace Head (Sipilä et al., 2016) and Arctic sites (Beck et al., 2021), the rate of the IA-MSA mechanism is also significantly lower. This indicates that the contribution of MSA to IA particle formation under atmospheric conditions is relatively limited, and a more efficient stabilizer for IA should be involved in the nucleating process, such as other iodine-containing components such as HIO₂ or iodine oxides (I₂O₄ and I₂O₅). Moreover, considering the complexity of the marine atmosphere, other non-iodine nucleation precursors, such as SA, NH₃, amines, etc., may also affect the nucleation process. This is particularly true with SA, because MSA and SA coexist in the air and both are formed during the oxidation of DMS in the marine atmosphere. Therefore, in future studies, the influence of the above factors on the nucleation mechanism of marine aerosols will also be considered.

4 Atmospheric significance and conclusion

The present work systematically investigates the joint nucleation mechanisms of two critical marine nucleation precursors, i.e., methanesulfonic acid (MSA) and iodic acid (IA), using the quantum chemical approach and Atmospheric Cluster Dynamics Code (ACDC). The results suggest that the self-nucleation rate of IA is much lower than that of the CLOUD experiment, indicating the importance of stabilizers for IA in the nucleation process. We find that MSA can stabilize IA cluster via both hydrogen and halogen bonds and thus promote IA cluster formation rate, especially in low-temperature environments with sparse IA and rich MSA. The corresponding IA-MSA nucleating mechanism can be described by two distinct pathways: (i) pure-IA cluster formation and (ii) IA-MSA cluster formation, indicating that IA and MSA can jointly nucleate. The IA-MSA nucleation is highly dependent on the distribution of MSA and IA in the marine atmosphere. However, IA-MSA nucleation rates are far from sufficient to explain the field observations, indicating that additional essential precursors need to be considered (e.g., H₂SO₄, HIO₂, NH₃, and amines). Nucleation among these components is likely to be synergistic, with only varying magnitudes of contribution. For example, both SA and MSA originate from the oxidation of DMS, so their coexistence in the atmosphere may synergistically promote the formation of IA clusters, which is worthy of future studies.

The current study provides molecular-level evidence that IA and MSA can jointly nucleate, and the IA-MSA joint nucleation is more efficient than the IA self-nucleation. More broadly, this finding helps to construct a more comprehensive marine multicomponent nucleation model.