the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Collisionsticking rates of acid–base clusters in the gas phase determined from atomistic simulation and a novel analytical interacting hardsphere model
Huan Yang
Ivo Neefjes
Valtteri Tikkanen
Jakub Kubečka
Theo Kurtén
Hanna Vehkamäki
Bernhard Reischl
Kinetics of collisionsticking processes between vapor molecules and clusters of lowvolatility compounds govern the initial steps of atmospheric new particle formation. Conventional noninteracting hardsphere models underestimate the collision rate by neglecting longrange attractive forces, and the commonly adopted assumption that every collision leads to the formation of a stable cluster (unit mass accommodation coefficient) is questionable for small clusters, especially at elevated temperatures. Here, we present a generally applicable analytical interacting hardsphere model for evaluating collision rates between molecules and clusters, accounting for longrange attractive forces. In the model, the collision cross section is calculated based on an effective molecule–cluster potential, derived using Hamaker's approach. Applied to collisions of sulfuric acid or dimethylamine with neutral bisulfate–dimethylammonium clusters composed of 1–32 dimers, our new model predicts collision rates 2–3 times higher than the noninteracting model for small clusters, while decaying asymptotically to the noninteracting limit as cluster size increases, in excellent agreement with a collisionratetheory atomistic molecular dynamics simulation approach. Additionally, we calculated sticking rates and mass accommodation coefficients (MACs) using atomistic molecular dynamics collision simulations. For sulfuric acid, a MAC ≈1 is observed for collisions with all cluster sizes at temperatures between 200 and 400 K. For dimethylamine, we find that MACs decrease with increasing temperature and decreasing cluster size. At low temperatures, the MAC ≈1 assumption is generally valid, but at elevated temperatures MACs can drop below 0.2 for small clusters.
 Article
(3831 KB)  Fulltext XML

Supplement
(965 KB)  BibTeX
 EndNote
In the first steps of atmospheric new particle formation (NPF), condensedphase clusters (Vehkamäki and Riipinen, 2012) form from supersaturated vapors (Guo et al., 2014; Kirkby et al., 2011; Kulmala et al., 2000, 2013; Zhang et al., 2012). Precise quantification of the cluster formation rates is important to accurately assess global aerosol budgets and their impact on climate (Hallquist et al., 2009; Ramanathan et al., 2001; Seinfeld et al., 2016) and health (Pope and Dockery, 2006; Pöschl, 2005). Cluster formation is conventionally approximated as a microscopic kinetic process where the size evolution of the newborn clusters is driven by competing events of monomer acquisition and loss. This results in a set of pseudofirstorder kinetic equations (Schenter et al., 1999) (written here for a onecomponent system for simplicity):
where N_{i} is the number concentration of clusters composed of i monomers (atoms or molecules), k_{i,1} the monomer–cluster collision rate coefficient, γ_{i} the evaporation rate coefficient, and α_{i,1} the mass accommodation coefficient which is a measure of the average probability for the monomer and cluster to “stick together” after a collision. The collision rate coefficient and mass accommodation coefficient can be combined into the sticking rate coefficient ${s}_{i,\mathrm{1}}={\mathit{\alpha}}_{i,\mathrm{1}}{k}_{i,\mathrm{1}}$. In principle, any cluster formation process can be accurately described by Eq. (1) and appropriate extensions (e.g., considering cluster dissociation or coagulation) if the associated rate coefficients are correctly determined. This work focuses particularly on the collision and sticking rate coefficients.
Precise determination of the collision and sticking rate coefficients is nontrivial, as they are sensitive to multiple variables including cluster size, structure, composition, system temperature, and vapor concentration. Conventional models, like the commonly used classical nucleation theory (CNT; Becker and Döring, 1935; Farkas, 1927), resort to crude assumptions to achieve generality. This inevitably sacrifices precision and leads to predictions agreeing with experiments only within a narrow range of saturation ratios and temperatures. In most conventional approaches, sticking rate coefficients are approximated by noninteracting hardsphere collision rate coefficients derived from kinetic gas theory. By neglecting intermolecular interactions, kinetic gas theory may not predict accurate collision rate coefficients, as intermolecular longrange attractive forces can enhance collisions, especially for polar or polarizable molecules common in atmospheric aerosol formation. Moreover, using collision rate coefficients to approximate sticking rate coefficients may be questionable, especially for small clusters, as the formation of a stable product cluster after a collision may be unsuccessful in a significant number of events; i.e., the mass accommodation coefficient can be lower than unity. The noninteracting hardsphere approximation thus potentially leads to significant errors in modeling the initial steps of cluster formation.
Molecular dynamics (MD) simulations allow us to study the time evolution of collision systems, where all intra and intermolecular interactions are described through force fields. Recently, we have developed several atomistic molecular dynamics simulation frameworks (Yang et al., 2018; Halonen et al., 2019; Neefjes et al., 2022) to determine collision probabilities from sampling trajectories of atoms or molecules, which predict enhanced collision rate coefficients in good agreement with experiments (Lehtipalo et al., 2016; Stolzenburg et al., 2020). However, as the adequate sampling of collision probabilities requires simulating a large number of binary collisions, these approaches are computationally expensive. More efficient ways need to be proposed to systematically calculate rate coefficients, in particular for collisions involving molecules and clusters of different sizes.
Rate coefficients can also be obtained through analytical models that consider the interactions between the collision partners. These models usually require a fraction of the computational cost of MD simulations but generally still rely on significant approximations. In the central field approach, the collision partners are approximated by point particles interacting through an isotropic attractive potential. Along this line, Gioumousis and Stevenson (1958) developed a model to calculate the secondary ion–molecule reaction rates in mass spectrometers. Su and Bowers (1973) derived a theory to describe ion–permanent dipole interactions and combined it with the central field approach to model ion–polar molecule collisions. Several other authors (Moran and Hamill, 1963; Su et al., 1978; Clary, 1985) have also theoretically studied collision rate coefficients between ions and molecules by applying the central field approach.
We have previously shown that the central field approach yields very similar results to those from sampling collision trajectories using MD, provided that the attractive interaction used in the central field model is fitted to the tail of the potential of mean force (PMF) between the collision partners (Neefjes et al., 2022). The largest collision partners considered in that study were, however, dimers. For collisions involving larger clusters, it becomes increasingly difficult to obtain accurate PMFs. Additionally, in the central field model, a collision is only considered to have occurred when the center of mass distance between the collision partners is infinitesimal. Such a collision criterion may underestimate cluster formation rates, as in realistic systems, cluster formation can already occur at larger centertocenter distances, where the extended structures of the collision partner begin to overlap.
In this study, we present two extensions to the central field approach to efficiently predict collision rate coefficients between monomers and clusters of arbitrary sizes.

We consider an interacting hardsphere model, in which the collision partners are still treated as point particles, interacting through an effective longranged isotropic attractive potential, but a collision is considered to have occurred if the center of mass distance between the collision partners is less than the sum of their hardsphere radii to account for their extended structures.

We integrate the monomer–monomer attractive potential following the approach of Hamaker (Hamaker, 1937) to obtain an effective monomer–cluster attractive potential for which a welldefined collision cross section (CCS) can be determined and the monomer–cluster collision rate coefficients can be calculated. Parameters of the monomer–monomer attractive potential are taken from the PMF between the collision partners, calculated from welltempered metadynamics simulations using atomistic force fields. To validate our analytical approach, we compared it to both atomistic MD collision simulations and MD simulations between two point particles interacting through the same effective potential.
Furthermore, we provide arguments from MD collision simulations for when a unit mass accommodation coefficient is valid and when this approximation fails. Sticking probabilities and sticking rate coefficients inferred from MD collision simulations can be used in Eq. (1) to model acid–baseinduced atmospheric NPF and compare it to experimental results.
As test systems, we considered collisions between the monomers of sulfuric acid (H_{2}SO_{4}) and dimethylamine ((CH_{3})_{2}NH) and clusters of dimethylammonium–bisulfate dimers ([${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n}, $n=\mathrm{1},\mathrm{2},\mathrm{4},\mathrm{8},\mathrm{16},\mathrm{32}$). Sulfuric acid is known to be a crucial component of continental NPF, while dimethylamine likely facilitates the formation of clusters through acid–base interactions (Kurtén et al., 2008; Sipilä et al., 2010). This paper is structured as follows: in Sect. 2.1, we present our analytical interacting hardsphere collision rate model and show how the CCSs can be calculated for monomer–monomer and monomer–cluster collisions using a simple van der Waals potential and an effective potential derived from Hamaker's approach. In Sect. 2.2, we describe the welltempered metadynamics calculations and the atomistic and point particle MD simulations. In Sect. 3, we validate the analytical model and present and discuss simulation results of the collision and sticking probabilities and rate coefficients for sulfuric acid–dimethylamine clusters. Section 4 summarizes and concludes the paper with an emphasis on its implications to atmospheric NPF modeling.
2.1 Analytical methods
2.1.1 Central field approach
In the central field model, both collision partners are represented by point particles. We consider the initial condition where one particle remains at rest while the other particle approaches from infinitely far away with some initial relative velocity v_{0}. The perpendicular distance between the initial velocity of the moving particle and the center of the particle at rest defines the impact parameter b. The geometry of an example system is illustrated in Fig. 1. This typical twodimensional twobody motion can be first reduced to a twodimensional singlebody orbital motion in a central field and then to a onedimensional singlebody motion along the particle centertocenter direction (Landau and Lifshitz, 1976), i.e.,
where ${U}_{\mathrm{eff}}\left(r\right)=U\left(r\right)+\mathit{\mu}{v}_{t}^{\mathrm{2}}/\mathrm{2}=U\left(r\right)+\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}{b}^{\mathrm{2}}/\mathrm{2}{r}^{\mathrm{2}}$ is an effective potential resulting from reducing the twodimensional orbital motion to the onedimensional motion along the centertocenter direction, U(r) the attractive potential between the point particles, r the centertocenter distance, v the instantaneous relative velocity, v_{r} the component of the instantaneous relative velocity along the centertocenter direction, v_{t} the component of the instantaneous relative velocity perpendicular to the centertocenter direction, and μ the reduced mass.
Reducing the twodimensional singlebody orbital motion to a onedimensional motion along the centertocenter direction introduces a centrifugal energy barrier (the repulsive part of the effective potential) between the point particles. In the central field model, a “collision” is defined by zero distance between the point particles. Crossing the centrifugal energy barrier in the centertocenter direction is, therefore, a necessary and sufficient condition for a collision. This condition is equivalent to the kinetic energy in the centertocenter direction $\mathrm{1}/\mathrm{2}\mathit{\mu}{v}_{\mathrm{r}}^{\mathrm{2}}$ being positive for all values of r>0; i.e., for a collision to happen, the following condition must hold for arbitrary r>0:
After substituting in the expression of the effective potential and rearranging the above equation, we obtain
where the function ω_{v}(r) is defined for convenience. As Eq. (4) must hold for arbitrary r>0, it is apparent that the minimum value of function ω_{v}(r) in interval r>0 sets an upper boundary for the impact parameter b, below which a collision will happen and above which a collision is not possible. This upper boundary defines the critical impact parameter b_{c,CF} in the central field model:
where ω_{v}(R_{m}) is the minimum of function ω_{v}(r) with R_{m} being the corresponding distance where the minimum is achieved. Physically, R_{m} is the closest distance that the collision partners can approach each other when the centrifugal barrier cannot be crossed.
Depending on the attractive potential U(r), Eq. (5) might, however, not provide a satisfactory definition for the critical impact parameter. Figure 2a shows the shape of ω_{v}(r) for general attractive potentials of the form $U\left(r\right)=A(r/{r}_{\mathrm{0}}{)}^{a}$. For this general form, we can distinguish three cases.

If $a<\mathrm{2}$, ω_{v}(r) is a convex function exhibiting a single minimum at $r={R}_{\mathrm{m}}={r}_{\mathrm{0}}\left[A\right(a+\mathrm{2})/\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}){]}^{\mathrm{1}/a}>\mathrm{0}$ and ${\mathit{\omega}}_{{v}_{\mathrm{0}},a<\mathrm{2}}\left({R}_{\mathrm{m}}\right)=a{r}_{\mathrm{0}}^{\mathrm{2}}\left[A\right(a+\mathrm{2})/\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}){]}^{\mathrm{2}/a}/(a+\mathrm{2})$.

If $a=\mathrm{2}$, ω_{v}(r) monotonically increases with a minimum at $r={R}_{\mathrm{m}}=\mathrm{0}$ and ${\mathit{\omega}}_{v,a=\mathrm{2}}\left({R}_{\mathrm{m}}\right)=\mathrm{2}A{r}_{\mathrm{0}}^{\mathrm{2}}/\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}$.

If $\mathrm{2}<a<\mathrm{0}$, ω_{v}(r) is a monotonically increasing function with a minimum at $r={R}_{\mathrm{m}}=\mathrm{0}$ and ${\mathit{\omega}}_{v,\mathrm{2}<a<\mathrm{0}}\left({R}_{\mathrm{m}}\right)=\mathrm{0}$.
These cases indicate that (1) for $\mathrm{2}<a<\mathrm{0}$, such as for the Coulombic potential, ${b}_{\mathrm{c},\mathrm{CF}}^{\mathrm{2}}={\mathit{\omega}}_{v,\mathrm{2}<a<\mathrm{0}}\left({R}_{\mathrm{m}}\right)=\mathrm{0}$ and the CCS is equally zero; the centrifugal barrier can never be crossed, and a collision is impossible regardless of the initial velocity and impact parameter. (2) For $a<\mathrm{2}$ and large initial velocities v_{0}, ${b}_{\mathrm{c},\mathrm{CF}}^{\mathrm{2}}={\mathit{\omega}}_{v,a<\mathrm{2}}\left({R}_{\mathrm{m}}\right)=a{r}_{\mathrm{0}}^{\mathrm{2}}\left[A\right(a+\mathrm{2})/\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}){]}^{\mathrm{2}/a}/(a+\mathrm{2})$, and the CCS likewise becomes vanishingly small. These issues arise because the collision partners are approximated as point particles in the central field approach; a collision is only considered to have occurred when the distance between the particles is zero. In realistic systems, the collision partners have an extended structure. As a result, the point particle approximation is expected to give a lower collision rate than when the extent of the system is considered.
2.1.2 Interacting hardsphere model
In the interacting hardsphere model, the aforementioned situations are accounted for by using a revised collision criterion; i.e., if the distance between the collision partners is at any point smaller than the sum of their hardsphere radii, R_{i}+R_{j}, a collision is considered to have occurred. It is apparent that when ${R}_{i}+{R}_{j}<{R}_{\mathrm{m}}$, this revised criterion will give an identical result as that in the central field model because in this case “reaching a distance smaller than R_{i}+R_{j}” is equivalent to “crossing the centrifugal barrier” for the point particles. The new criterion will make a difference only when R_{i}+R_{j}≥R_{m}, where even if the centrifugal barrier is not crossed, there is still a chance for the point particles to reach a distance smaller than R_{i}+R_{j} (see Fig. 2b). Using this revised criterion, we can define two cases for the critical impact parameter b_{c}:
The concept of determining the critical impact parameter by examining the aforementioned two cases has been discussed by other authors (Ouyang et al., 2012; Fuchs and Sutugin, 1965), but explicit expressions using the hardsphere radii, like Eq. (6), have not been presented in the existing literature to the best of our knowledge. Once the appropriate critical impact parameter is determined, the dynamical collision cross section (CCS) is given by the usual expression,
and the collision rate coefficient for interacting hard spheres is obtained as
where f(v_{0}) is the Maxwell–Boltzmann distribution. In the case of noninteracting hard spheres, the CCS is independent of the relative velocity, and the collision rate coefficient can be directly calculated by substituting $\mathrm{\Omega}\left(v\right)=\mathit{\pi}({R}_{i}+{R}_{j}{)}^{\mathrm{2}}$ into the above equation:
where k_{B} is the Boltzmann constant and T is the temperature.
To determine the critical impact parameter b_{c} in the interacting hardsphere model, we first need to find the value of R_{m} for a given interaction potential U(r), which will be done in the following sections.
Solution for the monomer–monomer potential
Vapors relevant to atmospheric cluster formation are mostly composed of polar and polarizable molecules interacting through a van der Waals potential (Leite et al., 2012) which scales with the separation distance as
where ϵ and σ are the characteristic energy and length of the attractive potential and the subscript “mm” denotes the monomer–monomer interaction. Substituting Eq. (10) into Eq. (4), the positive minimum of ω_{v}(r,U_{mm}) is found to be
located at
For a given pair of monomers, b_{c} can be readily evaluated from Eqs. (6) and (12) with respect to any given relative speed v_{0}. We obtained values for ϵ and σ from the potential of mean force (PMF) along the center of mass distance between the collision partners calculated from welltempered metadynamics calculations using an empirical atomistic force field, presented in Sect. 2.2 (results are included in the Supplement).
Solution for the monomer–cluster potential
Following the approach of Hamaker (Hamaker, 1937), an effective monomer–cluster potential U_{mc}(r) can be obtained by integrating the monomer–monomer potential over the volume of a spherical cluster, in which uniform monomer number density is assumed:
where ρ_{c} is the monomer number density, R_{c} the radius of the cluster, dV=ρ^{2}sin (ϕ)dθdϕdρ the volume element, V_{c} the total volume of the cluster, and n_{c}=ρ_{c}V_{c} the total number of monomers in the cluster. In Eq. (13), the effective monomer–cluster potential is essentially a sum over all individual contributions of the monomer–monomer potentials, which has two implications. First, the monomer–monomer potentials must be pairwise so that their individual contributions are additive. Hence, this approach is not suitable for nonpairwise interactions. Second, the approach can be conveniently extended to multicomponent clusters containing several types of monomers. By using Eq. (13) in a similar manner, the corresponding multicomponent monomer–cluster potential is $\mathrm{4}{\sum}_{i}{n}_{\mathrm{c},i}{\mathit{\u03f5}}_{i}{\mathit{\sigma}}_{i}^{\mathrm{6}}/{\left({r}^{\mathrm{2}}{R}_{\mathrm{c}}\right)}^{\mathrm{3}}$, where n_{c,i} is the number of monomers of type i in the cluster and ϵ_{i} and σ_{i} are the interaction parameters between a monomer of type i in the cluster and the free monomer colliding with the cluster. Note that further integrating U_{mc}(r) over the volume of another cluster recovers the wellknown result of Hamaker for a van der Waals attractive potential between two spherical objects, which can be invoked to treat cluster–cluster collision rates. Integrated cluster–cluster potentials have previously been used to model collision rates, but the models usually resort to fitting against experimental data (Fuchs and Sutugin, 1965; Sceats, 1989). Replacing U(r) in Eq. (4) with U_{mc}(r), we find ω_{v}(r,U_{mc}) to be a convex function with a positive minimum (note that we require r>R_{c}). By taking the derivative of ω_{v}(r,U_{mc}) with respect to r, the minimum is found to be located at r=R_{m}(U_{mc}) and satisfies the following:
with ${a}_{\mathrm{0}}={R}_{\mathrm{c}}^{\mathrm{2}}({R}_{\mathrm{c}}^{\mathrm{6}}{l}_{\mathrm{c}}^{\mathrm{6}})$, ${a}_{\mathrm{1}}=\mathrm{2}{l}_{\mathrm{c}}^{\mathrm{6}}\mathrm{4}{R}_{\mathrm{c}}^{\mathrm{6}}$, ${a}_{\mathrm{2}}=\mathrm{6}{R}_{\mathrm{c}}^{\mathrm{4}}$, ${a}_{\mathrm{3}}=\mathrm{4}{R}_{\mathrm{c}}^{\mathrm{2}}$, and a_{4}=1, where ${l}_{\mathrm{c}}\mathrm{\equiv}[\mathrm{8}{n}_{\mathrm{c}}\mathit{\u03f5}{\mathit{\sigma}}^{\mathrm{6}}/(\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}){]}^{\mathrm{1}/\mathrm{6}}$ is a cluster size and potentialdependent characteristic length. Equation (14) is a fourthorder (quartic) equation after writing it in terms of ${R}_{\mathrm{m}}^{\mathrm{2}}$, and it can be shown that only one real root exists for R_{m}>R_{c} (see Supplement). The solution is readily available:
where $q=\mathrm{2}{l}_{\mathrm{c}}^{\mathrm{6}}$, $M=\sqrt{(N+{\mathrm{\Delta}}_{\mathrm{0}}/N)/\mathrm{3}}/\mathrm{2}$, and $N=\sqrt[\mathrm{3}]{({\mathrm{\Delta}}_{\mathrm{1}}+\sqrt{{\mathrm{\Delta}}_{\mathrm{1}}^{\mathrm{2}}\mathrm{4}{\mathrm{\Delta}}_{\mathrm{0}}^{\mathrm{3}}})/\mathrm{2}}$ with ${\mathrm{\Delta}}_{\mathrm{0}}=\mathrm{36}{R}_{\mathrm{c}}^{\mathrm{2}}{l}_{\mathrm{c}}^{\mathrm{6}}$ and ${\mathrm{\Delta}}_{\mathrm{1}}=\mathrm{108}{l}_{\mathrm{c}}^{\mathrm{12}}$. With Eqs. (6) and (15), we can predict b_{c} for any monomer and cluster pair, as long as the monomer–monomer interaction parameters ϵ and σ are provided.
2.2 Computational methods
Computer simulations were used to

compute the PMF between monomer pairs to obtain the ϵ and σ parameter values required for the analytical collision rate coefficient model;

validate the analytical model by comparison with atomistic and point particle molecular dynamics (MD) collision simulations;

calculate predictions for the sticking probability, sticking rate coefficient, and mass accommodation coefficient by analyzing the atomistic MD collision simulations over a range of relative velocities and impact parameters.
2.2.1 Systems and force field
We considered the atmospherically relevant sulfuric acid (H_{2}SO_{4}) and dimethylamine ((CH_{3})_{2}NH) molecules, as well as clusters consisting of 1, 2, 4, 8, 16, or 32 bisulfate–dimethylammonium dimers ([${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]). Clusters composed of n [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] dimers were obtained by sintering and equilibrating two smaller clusters with $n/\mathrm{2}$ dimers.
In the PMF and atomistic MD simulations, these systems are described by a force field fitted according to the OPLS (optimized potentials for liquid simulations) allatom procedure (Jorgensen et al., 1996). In the OPLS force field, the intramolecular interactions consist of harmonic bond potentials between covalently bonded atoms, harmonic angle potentials between atoms separated by two covalent bonds, and dihedral angle potentials between atoms separated by three covalent bonds,
where ${k}_{i}^{\mathrm{b}}$, r_{i}, and ${r}_{i}^{\mathrm{0}}$ are the force constant, instantaneous, and equilibrium length of bond i; ${k}_{j}^{\mathit{\theta}}$, θ_{j}, and ${\mathit{\theta}}_{j}^{\mathrm{0}}$ are the force constant, instantaneous, and equilibrium value of angle j; and V_{n}, ${\mathit{\varphi}}_{n}^{k}$, and ϕ^{k} are the Fourier coefficients, phase angles, and instantaneous value of the dihedral angle k.
The intermolecular interactions, as well as intramolecular interactions between atoms separated by more than three covalent bonds, are described by LennardJones potentials between atoms i and j separated by a distance r_{ij}, with distance and energy parameters σ_{ij} and ϵ_{ij}, as well as Coulomb interactions between the atoms' partial charges q_{i} and q_{j}:
where ϵ_{0} is the vacuum permittivity.
The OPLS force field parameters used in this study were obtained from Loukonen et al. (2010). We note that in the original OPLS force field, LennardJones and Coulomb interactions between atoms separated by three covalent bonds (“1–4 interactions”) are scaled by a factor 0.5. Loukonen et al. (2010) set this scaling factor to zero when fitting the force field parameters. For consistency, we have also set these interactions to zero in our simulations.
We have validated these force field parameters for MD collision simulations in our previous work (Halonen et al., 2019) by comparing predicted molecular structures, binding energies, and vibrational spectra of sulfuric acid molecules with ab initio results.
2.2.2 Potential of mean force calculations
To determine the thermally averaged interaction potential between the collision partners, we calculated the PMF as a function of their center of mass distance from welltempered metadynamics simulations (Barducci et al., 2008), using the PLUMED plugin (Tribello et al., 2014) for LAMMPS (Plimpton, 1995). Calculations were performed for five combinations of monomer pairs at 300 K: (CH_{3})_{2}NH–(CH_{3})_{2}NH, H_{2}SO_{4}–H_{2}SO_{4}, H_{2}SO_{4}–(CH_{3})_{2}NH, H_{2}SO_{4}–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1}, and (CH_{3})_{2}NH–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1}. We used the velocity Verlet integrator with a time step of 1 fs, where the LennardJones interactions were cut off at 14 Å, while electrostatic interactions were evaluated with a cutoff at 40 Å. A total of 40 random walkers were employed to deposit Gaussians with a width of 0.1 Å and initial height of 2k_{B}T (${k}_{\mathrm{B}}T/\mathrm{10}$ for the (CH_{3})_{2}NH–(CH_{3})_{2}NH system) every 500 time steps along the collective variable; a harmonic wall was used to restrict the collective variable to values below 35 Å. A bias factor of 20 (5 for the (CH_{3})_{2}NH–(CH_{3})_{2}NH system) was chosen, and a stochastic velocity rescaling thermostat with a time constant of 0.1 ps was used to maintain a constant temperature. Note that we used the PMF at 300 K to approximate cases at 200 and 400 K in deriving corresponding collision rate coefficients. Although there are slight differences in the PMFs at these temperatures, their influence on the analytical collision rate coefficients was found to be negligible. Results are summarized in Fig. S1 and Table S1 of the Supplement. The attractive interaction between (CH_{3})_{2}NH and (CH_{3})_{2}NH was too small to be distinguished from thermal noise in our calculations. It was therefore directly set to zero and is not shown in Fig. S1. The monomer–monomer interactions were assumed to take the form of a 12–6 LennardJones potential. The depth of the potential well (or binding energy) ϵ was directly taken as the minimum energy of the PMF curve, and σ was determined from $\mathit{\sigma}={r}_{\mathit{\u03f5}}/\sqrt[\mathrm{6}]{\mathrm{2}}$, where r_{ϵ} is the corresponding distance for the minimum energy point on the PMF curve. The binding energy for the H_{2}SO_{4}–H_{2}SO_{4} system was found to be −0.29 eV, which is in excellent agreement with previous ab initio results (Temelso et al., 2012) of −0.3 eV obtained from Boltzmann averaging over four energy minimum structure dimers at 298.15 K and is in good agreement with recent calculations (Elm et al., 2016; Myllys et al., 2017) of −0.23 and −0.26 eV at a higher level of theory. As the outcome of the collisions is already decided at a significant separation between the collision partners for the systems considered here, the role of proton transfer after a successful collision and other quantum effects can be ignored, and the classical force field approach is expected to capture the most important features of the collision.
We note that, though in the current paper the monomer–monomer interacting parameters were obtained from the PMF calculated by atomistic simulations, they could in principle be taken directly from literature values if available.
2.2.3 Atomistic molecular dynamics collision simulations
Atomistic simulations were performed with the LAMMPS code (Plimpton, 1995) for collisions of either H_{2}SO_{4} or (CH_{3})_{2}NH with [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n}, where $n=\mathrm{1},\mathrm{2},\mathrm{4},\mathrm{8},\mathrm{16}$, or 32 at temperatures of 200, 300, and 400 K. Equations of motion were integrated using a velocity Verlet algorithm with a time step of 1 fs, where LennardJones and electrostatic interactions were evaluated with a cutoff at 280 Å. The collision partners are originally placed 550 Å apart along the x axis, well beyond the cutoff of the intermolecular potentials. At these positions, the collision partners were equilibrated for 50 ps using separate Langevin thermostats for each collision partner with a damping factor of 0.1 ps. The atomic velocities were drawn from the Maxwell–Boltzmann distribution at the target temperature with the center of mass motion of each collision partner and angular momentum of the combined system removed. After equilibration, the distance between the now orientationally randomized collision partners was decreased to 150 Å along the x axis, bringing them within range of the longrange intermolecular potentials. Both collision partners were then given a velocity along the x direction of ${v}_{x}=\pm {v}_{\mathrm{0}}/\mathrm{2}$ towards each other, where v_{0} is the initial relative velocity. For each system, the range of relative velocities started at 50 m s^{−1} and increased in steps of 50 m s^{−1}. The highest relative velocity was determined so that at least 99 % of the Maxwell–Boltzmann distribution was sampled. We sampled impact parameters starting from 0 up to 36 Å, in steps of 1 Å along the z axis. At initial relative velocities of 50 m s^{−1}, there is still some probability for a collision even beyond an impact parameter of 36 Å. These probabilities, however, do not contribute significantly to the collision rate coefficient.
For the collision systems with clusters containing 1, 2, 4, 8, 16, and 32 [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] dimers, we performed N=1000, 1000, 800, 400, 200, and 100 independent trajectory simulations, respectively, for each combination of relative velocity and impact parameter. These result in statistically significant estimations of the collision probability. The resulting trajectories were analyzed to determine the collision probability: if the collision partners’ center of mass distance was smaller than R_{i}+R_{j} (determined once again from the mass and bulk density by assuming spherical shapes) in at least one time frame, then a collision event was identified. This definition is consistent with that in the analytical model, and it leads to the collision probability ${P}_{\mathrm{c}}(b,{v}_{\mathrm{0}})={N}_{\mathrm{c}}/N$, where N_{c} is the numbers of identified collision events. The definition of the sticking probability is similar to that of the collision probability, i.e., ${P}_{\mathrm{s}}(b,{v}_{\mathrm{0}})={N}_{\mathrm{s}}/N$, where the criterion for determining N_{s} is discussed in detail in Sect. 3.3.1.
2.2.4 Point particle molecular dynamics collision simulations
In the point particle molecular dynamics (PP MD) simulations, collision partners were treated as point masses interacting through effective pair potentials. The purpose of performing such simulations was to verify the critical impact parameters derived from the analytical interacting hardsphere model. The collision simulation setup in this case is similar to that in the atomistic MD simulation. First, the point particles were placed 200 Å apart on the x axis, with a certain impact parameter along the y axis, and given an initial relative velocity along the x direction. Second, the equations of motion of the point particles were solved with the velocity Verlet algorithm. A dynamic time step interval was used to ensure that the displacement of the point particles at each time step cannot exceed 0.1 Å. A collision simulation was stopped when either of the following two situations happened: (1) the distance between the point particle pair is less that the sum of their hardsphere radii (given in the analytical interacting hardsphere model), or (2) the distance between the point particle pair is larger than the starting distance. Situation 1 corresponds to a collision event, while situation 2 indicates that a collision will not happen in that simulation. The critical impact parameter for a given initial relative speed was found by gradually increasing the impact parameter (starting from the summation of hardsphere radii given in the interacting hardsphere model) until the transition from situation 1 to situation 2 was observed. This procedure was used to simulate two example cases where point particle pairs interact through either a Coulomb or a LennardJones potential. Simulated critical impact parameters as a function of initial relative speeds and their comparisons with those predicted from the analytical interacting hardsphere model are discussed in Sect. 3.1.
2.2.5 Collision rate, sticking rate, and mass accommodation coefficients
The MD collision rate coefficient k_{MD} was calculated by numerically integrating the collision probability as a function of impact parameter and relative initial velocity obtained from the MD trajectory simulations over the impact parameter and relative speed distribution at the target temperature:
Note that Eq. (18) reduces to Eq. (8) by letting ${P}_{\mathrm{c}}(b,{v}_{\mathrm{0}})=\mathrm{1}$ for b<b_{c} and ${P}_{\mathrm{c}}(b,{v}_{\mathrm{0}})=\mathrm{0}$ for b>b_{c}. Similarly, the sticking rate coefficient is calculated through
The mass accommodation coefficient is defined as
and it characterizes the average probability of sticking after collisions at the specified temperature.
3.1 Validation of the analytical interacting hardsphere model
To validate our analytical interacting hardsphere model, we first compared the analytical critical impact parameter b_{c} predicted by Eq. (6) to molecular dynamics simulations of two point particles interacting through the same effective pair potential (i.e., PP MD). Two cases were considered: (1) collisions between the $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$ ion and ${\mathrm{HSO}}_{\mathrm{4}}^{}$ ion under the attractive Coulomb potential $U={e}^{\mathrm{2}}/\left(\mathrm{4}\mathit{\pi}{\mathit{\u03f5}}_{\mathrm{0}}r\right)$, where e is the elementary charge and ϵ_{0} the vacuum permittivity (Fig. 3a), and (2) collisions between the neutral sulfuric acid H_{2}SO_{4} monomer and the neutral acid–base [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1} dimer under the van der Waals potential $U=\mathrm{4}\mathit{\u03f5}(\mathit{\sigma}/r{)}^{\mathrm{6}}$, with parameters obtained from the potential of mean force (PMF) along the center of mass (COM) distance between the molecules (Fig. 3b). In the PP MD simulations, b_{c} was found for a specific relative speed by gradually increasing b until the distance r between the particles is never smaller than the sum of the hardsphere radii. In the analytical model, R_{m}=0 for the attractive Coulomb potential and ${R}_{\mathrm{m}}=\mathit{\sigma}{\left[\mathrm{16}\mathit{\u03f5}/\left(\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}\right)\right]}^{\mathrm{1}/\mathrm{6}}$ for the van der Waals potential.
Excellent agreement is observed in both cases (note that the dependent variable b_{c} is set as the x axis in Fig. 3). For the Coulomb potential, R_{i}+R_{j} is always larger than R_{m}, and ${b}_{\mathrm{c}}={\mathit{\omega}}_{v}({R}_{i}+{R}_{j}{)}^{\mathrm{1}/\mathrm{2}}$ in the analytical model. The Coulomb potential case discussed here is only used to verify Eq. (6) but should not be implemented further in Eqs. (7) and (8) to calculate collision rate coefficients, as the Coulomb potential is nonnegligible at distances comparable to the mean free path of the colliding ions at atmospheric pressures and hence violates the assumption of free molecular regime. For the van der Waals potential, ${R}_{\mathrm{m}}>{R}_{i}+{R}_{j}$ at low relative speeds, but at relative velocities over ∼ 800 m s^{−1} the situation is reversed and ${R}_{\mathrm{m}}<{R}_{i}+{R}_{j}$. For velocities above ∼ 800 m s^{−1}, the analytical b_{c} transitions from ${\mathit{\omega}}_{v}({R}_{\mathrm{m}}{)}^{\mathrm{1}/\mathrm{2}}$ to ${\mathit{\omega}}_{v}({R}_{i}+{R}_{j}{)}^{\mathrm{1}/\mathrm{2}}$. It is important to note that for the H_{2}SO_{4}–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1} system, relative speeds over ∼ 800 m s^{−1} account for less than 1 % of the Boltzmann relative speed distribution. In this case, ${b}_{\mathrm{c}}={\mathit{\omega}}_{v}({R}_{\mathrm{m}}{)}^{\mathrm{1}/\mathrm{2}}$ is a good approximation. However, as ${R}_{\mathrm{m}}\propto \mathit{\sigma}(\mathit{\u03f5}/\mathit{\mu}{v}_{\mathrm{0}}^{\mathrm{2}}{)}^{\mathrm{1}/\mathrm{6}}$, smaller values of ϵ and σ or larger values of ω and v_{0} can cause this transition to occur at lower relative speeds and increase the importance of the ${b}_{\mathrm{c}}={\mathit{\omega}}_{v}({R}_{i}+{R}_{j}{)}^{\mathrm{1}/\mathrm{2}}$ case in the calculation of collision rate coefficients. In the Supplement, we show examples of collisions of carbon dioxide and water molecules, where considering ω_{v}(R_{i}+R_{j}) is important to obtain reasonable collision cross sections. These systems will be studied in further detail in a future publication.
In Fig. 4, we compare the critical impact parameters obtained from the analytical model (Eq. 6) to the collision probability P_{c}(b,v_{0}) obtained from atomistic collision MD simulations, sampling the relevant range of impact parameters and relative velocities. We examine both “monomer–monomer” (H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1}) and “monomer–cluster” (H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{16}) collisions. For the former, the effective interaction in the analytical model is obtained from the PMF as a function of the monomer–monomer COM distance (see Fig. S1 and Table S1 in the Supplement). For the latter, the effective interaction is obtained from Hamaker's approach (Eq. 13). Note that the atomistic force field describing the molecules and ions in the collision MD simulations is the same as in the PMF calculation, from which we derive the effective interaction parameters used to calculate the analytical critical impact parameters.
In general, P_{c}(b,v_{0}) exhibits a clear boundary between values of b and v_{0} where collisions are highly probable (white region) and values where the probability is essentially 0 (black region). The location of this boundary depends on the intermolecular forces between the collision partners. An exception to this clear boundary can be seen at low relative speeds. Here, ${P}_{\mathrm{c}}(b,{v}_{\mathrm{0}})<\mathrm{1}$ even for small values of b, while at the same time, P_{c}(b,v_{0}) remains nonzero even for large values of b. This issue arises from the asymmetric structure of the molecules: due to the low relative speed, small but periodic repulsive interactions at certain collision angles can result in the molecules rebounding from each other before the orbit crosses the sum of the hardsphere radii (Halonen et al., 2019). These low relative speeds, however, only account for a very small percentage of the Boltzmann speed distribution.
For both systems, the analytical b_{c} values obtained from ω_{v}(R_{m}) are in excellent agreement with the collision probabilities obtained from atomistic collision MD simulations over the relevant range of relative velocities at 300 K. The transition to ${R}_{\mathrm{m}}>{R}_{i}+{R}_{j}$ occurs for relative velocities around 800 ms^{−1} for H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1} and around 2300 m s^{−1} for H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{16}. Once again, this indicates that the ${b}_{\mathrm{c}}={\mathit{\omega}}_{v}({R}_{i}+{R}_{j}{)}^{\mathrm{1}/\mathrm{2}}$ case is not relevant for these systems at the studied temperature of 300 K. Most importantly, the excellent agreement between analytical b_{c} values and collision probabilities from atomistic collision MD simulations for the “monomer–cluster” collision demonstrates the validity of Hamaker's approach in modeling the effective monomer–cluster interactions of the studied systems.
3.2 Collision rate enhancement factors predicted by the analytical interacting hardsphere model
Having validated the analytical interacting hardsphere model by comparison with MD simulations, we now use the model to predict enhancement factors η over the noninteracting hardsphere model. We define $\mathit{\eta}={k}_{\mathrm{IHS}}/{k}_{\mathrm{HS}}$ as the ratio between the collision rate coefficient obtained from the interacting and noninteracting hardsphere model. Figure 5 shows η for H_{2}SO_{4}–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} collisions and (CH_{3})_{2}NH–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} collisions at 300 K. Both monomer types have similar η: while the H_{2}SO_{4} monomer has stronger intermolecular interactions, the (CH_{3})_{2}NH monomer has larger thermal mean speeds. For both cases, the monomer–cluster collisions gradually approach noninteracting hardsphere behavior as the size of the cluster increases, i.e., η→1 as n→∞. At smaller size ranges the noninteracting hardsphere model can underestimate the collision rate coefficient by a factor of 2 to 3, potentially introducing nonnegligible systematic errors in atmospheric NPF models. Hence, the correction due to longrange attractive forces may need to be included in NPF models.
3.3 Analysis of molecular dynamics collision simulations
3.3.1 Timescale analysis and issues related to classical force field models
In this work, we use classical, atomistic molecular dynamics simulations to determine parameters used in the analytical model and to validate the analytical model by comparing them to collision MD simulations. It is therefore important to recall the limitations of classical force fields regarding the physicochemical processes occurring in gasphase collisions and the timescales on which these different processes occur.
Whether or not a collision will occur is typically determined at a relatively large intermolecular distance, where the collision partners can accurately be modeled using classical force field descriptions in a molecular dynamics simulation. Modeling the probability of the collision partners forming a stable product using MD simulations is, however, significantly more complex, as the “sticking” probabilities are influenced by the possibility of bond formation between the collision partners, which is not captured by nonreactive classical force fields. The acid–base clusters considered in this study can form hydrogen bonds upon collision, which is captured by the classical force field employed in the atomistic simulation; however, proton transfer between acids and bases cannot happen.
After a collision, the newly formed cluster will generally have some excess energy due to these bond formations. In the atmosphere, this energy will be removed from the cluster through collisions with surrounding gas molecules. In the collision MD simulations, however, no additional carrier gas is modeled. As such, the excess energy after a collision cannot be released.
Lastly, it is likely that the cluster dissociation process is unphysically enhanced due to the fact that in the classical atomistic model employed, each vibrational mode possesses ${k}_{\mathrm{B}}T/\mathrm{2}$ energy on average, while in a quantummechanical description of the same system, some highfrequency degrees of freedom would remain “frozen” at atmospherically relevant temperatures.
The mentioned limitations of MD simulations in modeling the sticking probability all enhance the dissociation of the formed cluster. As such, the sticking probability obtained from MD simulations can be seen as a limiting value. If the collision partners are able to stick together for a certain time after collision in the MD simulations, then it is likely that the same collision partners would be able to stick together for the same amount of time in the atmosphere. MD simulations are, however, limited in the timescale that they can model due to the computational cost. We, therefore, discuss some characteristic timescales in atmospheric cluster formation after collision.
We first differentiate the two types of monomer–cluster interactions influencing the cluster formation process: (1) the vibration and diffusion of a condensable vapor monomer on the cluster surface immediately after a collision, which, if successful, drives the cluster formation while storing excess energy due to bond formation (hydrogen bonds in this study), and (2) collisions between background carrier gas and the cluster, dissipating the excess energy and equilibrating/thermalizing the nascent cluster. We consider the following scenario: after a condensable vapor–cluster collision, the vapor monomer experiences a few rounds of unsteady vibrations on the cluster surface before forming stable hydrogen bonds (i.e., sticking) or rebounding from the cluster surface. In the former case, a stable cluster with excess thermal energy, distributed over its degrees of freedom, is produced. At this point, the captured monomer can only escape from the cluster through thermal fluctuation (i.e., evaporation), as the monomer has become thermally indistinguishable from other molecules in the cluster. At typical atmospheric conditions, the timescale t_{s} for the formation of stable bonds only requires a few vibrations of the impinging monomer on the cluster surface. Hence, t_{s} is commonly smaller than the mean free time t_{c} for the carrier gas–cluster collision and the time t_{eq} required to fully thermalize the cluster. Moreover, t_{s}, t_{c}, and t_{eq} can all be assumed to be orders of magnitude smaller than the time t_{ev} for a monomer to evaporate (see the schematic diagram in Fig. 6).
For example, for H_{2}SO_{4}–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] and (CH_{3})_{2}NH–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] surrounded by air at standard atmospheric pressure and 300 K, ${t}_{\mathrm{s}}\sim \mathrm{2}{L}_{\mathrm{c}}/\overline{v}\sim \mathrm{10}$ ps, estimated based on its molecular neighbor separation distance L_{c} and the mean thermal speed $\overline{v}$. Based on the molecular collision frequency, t_{c}∼10^{2} ps. As thermalization requires tens to hundreds of collisions (Yang et al., 2022), t_{eq}∼10^{3} ps. For these systems, by far the shortest evaporation timescale t_{ev}∼10^{6} ps is observed for the evaporation of a (CH_{3})_{2}NH from the (CH_{3})_{2}NH[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] cluster. The timescales of all other possible evaporation processes are significantly longer, ${t}_{\mathrm{ev}}\sim {\mathrm{10}}^{\mathrm{17}}{\mathrm{10}}^{\mathrm{22}}$ ps (Ortega et al., 2012). Therefore, for the systems considered here, even the shortest t_{ev} is significantly larger than t_{s}, t_{c}, and t_{eq}.
The timescale analysis strongly indicates that a stable cluster is very likely to form if the impinging condensable vapor monomer can survive a very short window of unsteady vibration right after a collision.
3.3.2 Collision and sticking criteria
To calculate collision and sticking probabilities in atomistic molecular dynamics simulations we need simple, but robust, criteria for when two collision partners have collided and when they are sticking together, as the large number of individual trajectories makes a manual inspection unfeasible. Collision and sticking criteria can be based on geometric and/or energetic properties of the system (Yang et al., 2018; Halonen et al., 2019; Neefjes et al., 2022). Here, we considered a criterion based on the center of mass distance between the collision partners, as this does not require knowledge of the detailed bonding between monomers and clusters of different sizes and can be directly related to the analytical interacting hardsphere model.
In analyzing the MD collision simulations, we broadly encounter four different kinds of trajectories: flyby, rebounding collision, collision with metastable product, and collision with stable product. An example of each of these types of trajectories is shown in Fig. 7. A flyby is characterized by a clearly parabolic curve when plotting the center of mass distance between the collision partners as a function of time. The curves of all collision trajectories show a bend towards a faster decrease in center of mass distance over time due to attractive interactions between the collision partners. In rebounding collisions, the trajectory sharply reverses direction after reaching its lowest intermolecular distance. In the two other collision types, a cluster is formed and the trajectory oscillates around some intermolecular distance. The collision partners can either dissociate after some time (metastable) or stay together longer than the assumed timescale for cluster equilibration (stable). As the impact parameter decreases, it becomes more difficult to distinguish between a flyby and a collision.
A distance criterion that is too strict runs the risk of miscounting collisions with relatively large center of mass distances arising from certain relative orientations as flybys. Conversely, a distance criterion that is too loose results in certain close flybys being counted as rebounding collisions. Generally, a good distance criterion would lie just above the largest amplitudes of the oscillating sticking trajectories in Fig. 7 and would need to be found for each combination of T, v_{0}, and b, which is not possible due to the large number of combinations that are analyzed. We found the sum of hardsphere radii of the studied collision partners presented a surprisingly good distance criterion for the studied systems over a wide range of impact parameters and relative velocities. A comparison between the hardsphere radii distance criterion and a manually optimized distance criterion yielded differences in collision probabilities of at most 1 %.
Based on our analyses, a collision is considered to have occurred if the center of mass distance between the collision partners is less than the sum of their hardsphere radii in at least one time frame of the trajectory. Based on our timescale analyses in Sect. 3.3.1, a formed cluster is considered “sticking” if the center of mass distance between the collision partners is less than the sum of hardsphere radii at any time between 25 ps after collision and the end of the simulation.
3.4 Collision and sticking rate coefficients from molecular dynamics simulations
Based on the above analysis, we monitored the center of mass distance between the impinging monomer and cluster during the atomistic MD simulations. Collision probabilities and sticking probabilities 25 ps or 50 ps after collision are shown in Fig. 8 for two different cases: (A) collisions between H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1} and (B) collisions between (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1}. The figure also shows the 1000 underlying trajectories at the relative velocity and impact parameter closest to the mode of the Maxwell–Boltzmann distribution and the sum of hardsphere radii, respectively, for both systems. There are several implications of the results in Fig. 8. First, the collision and sticking probabilities are almost identical for the system of H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1}, suggesting that most collisions lead to the formation of stable product clusters. For the system of (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1}, the sticking probability is visibly lower than the collision probability, indicating that a significant portion of formed clusters dissociate in the first 25 ps after collision. This can largely be attributed to the fact that the hydrogen bond between the (CH_{3})_{2}NH monomer and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} cluster is weaker than that between the H_{2}SO_{4} monomer and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} cluster. Second, even though some dissociations can still be observed in the underlying trajectories of the system of (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{1} later than 25 ps after collision, the sticking probability at 25 ps is almost identical to that at 50 ps for both systems, indicating that increasing t_{s} beyond 25 ps has little effect on the sticking probability; 25 ps is sufficient for the formation of stable product clusters, which is consistent with our characteristic timescale analysis.
The sticking rate coefficients and rate enhancement factors for different clusters sizes and different temperatures are shown in Fig. 9a for H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} collisions and Fig. 9b for (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} collisions, where circles and squares correspond to results obtained at 25 and 50 ps, respectively. The sticking rate enhancement factor is an ad hoc quantity defined as the ratio between the sticking rate coefficient and the hardsphere collision rate coefficient. For comparison, the collision rate coefficients from both the analytical model (solid line) and MD simulations (stars) are also shown in Fig. 9. For collisions between H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} (Fig. 9a), the collision rate coefficients and the sticking rate coefficients obtained at the two different times are all identical, as the hydrogen bonds are strong enough to “capture” all impinging monomers colliding with the cluster. For collisions between (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n} (Fig. 9b), several interesting trends can be observed. At a low system temperature of 200 K, the collision rate coefficients are almost identical to the sticking rate coefficients except for a few data points at small cluster sizes. At this temperature, the binding energy of the hydrogen bond is large enough compared to the low kinetic energy of the impinging monomer to prevent rebounds. However, the deviation between the collision and sticking rate coefficients becomes visible with increasing system temperature, which leads to higher kinetic energies of the impinging monomer, resulting in more rebounding. At a system temperature of 300 K, this deviation is expected to be solely a consequence of the unsuccessful vibrational coupling events (i.e., rebounding), whereas at T=400 K, the characteristic evaporation time t_{ev} of the nascent product cluster has decreased to a value comparable to t_{s}. At 400 K, the deviation between collision and sticking rate coefficients is caused by the coupled effects of the unsuccessful vibrational coupling events and the evaporation of the nascent product cluster. Consequently, visible differences between sticking rate coefficients obtained at 25 and 50 ps appear at 400 K. Moreover, the deviation between sticking and collision rate coefficients decreases with increasing cluster size: larger sizes favor the formation of stable clusters. This is not surprising as larger clusters offer more potential sites for hydrogen bonds and more molecular neighbors to enhance mutual attractions. They also serve as more efficient sinks for dissipating the latent heat released during hydrogen bond formation upon collision (Yang et al., 2019, 2022).
It is important to remember the limitations of MD simulations using a classical force field description and the timescale analysis performed in Sect. 3.3.1 when predicting sticking rate coefficients. Overall, we expect our simulation approach to overestimate the probability of dissociation. Based on our simulations, it is therefore very likely that the H_{2}SO_{4}–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] system indeed has a mass accommodation coefficient near unity. For the (CH_{3})_{2}NH–[${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] system, our simulations suggest a mass accommodation coefficient significantly less than unity at room temperature and above. Based on the weaker interactions between (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$] this is plausible, but the limitations of our simulation approach do not permit us to make a definitive statement. Ideally, collision would be studied using ab initio or at least reactive MD simulations in the presence of carrier gas, and for much longer timescales, to properly model all postcollision processes. Such calculations are, however, beyond the scope of current computational possibilities.
Precise rate coefficients are essential in modeling the gasphase collisions which are the initial step in atmospheric new particle formation (NPF). The noninteracting hardsphere model commonly adopted for neutral molecules and clusters may significantly underestimate collision rate coefficients due to neglecting longrange attractive forces. On the other hand, the assumption of a unit mass accommodation coefficient can overestimate the number of stable clusters formed, in particular for small clusters at elevated temperatures. It is important to note, however, that these two errors will not necessarily cancel each other out.
Simulation approaches based on atomistic modeling of the collision partners have recently been developed and shown to give reasonable estimates for collision rates and enhancement factors for collisions of individual molecules or ions. To systematically calculate rate coefficients for collisions between molecules and clusters of atmospherically relevant acid–base systems, we have developed a new analytical interacting hardsphere model. The molecule–cluster interactions are obtained from Hamaker's approach by integrating monomer–monomer interactions over the volume of the cluster. Here, the underlying monomer–monomer interacting parameters were obtained from fitting LennardJones potential to the monomer–monomer potential of mean force calculated from atomistic simulations, but we note that the monomer–monomer interacting parameters could also be obtained by other methods or taken directly from literature values if available. The critical impact parameter is determined using the sum of hardsphere radii or the minimum distance between the pointlike collision partners which cannot be crossed due to the centrifugal barrier in the central field model. The accuracy of the analytical model was validated against molecular dynamics collision simulations using the full atomistic model of the collision partners. The analytical model has an accuracy comparable to atomistic molecular dynamics simulations but can be applied efficiently for the systematic calculation of molecule–cluster collision rates required in atmospheric NPF models. The same approach can also be applied to calculate effective cluster–cluster interactions.
For collisions of sulfuric acid (H_{2}SO_{4}) or dimethylamine ((CH_{3})_{2}NH) molecules with clusters consisting of 132 bisulfate–dimethylammonium ([${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]) dimers, we find that for cluster radii smaller than 2 nm, the noninteracting hardsphere model underestimates the collision rate coefficient by a factor of 2 to 3 due to the neglect of the longrange attractive forces. The enhancement factor drops below 1.2 for cluster radii larger than 25 nm. This deviation from the noninteracting hardsphere model can be large enough to introduce systematic errors in atmospheric NPF models. For all the systems considered in this study, the critical impact parameter was well described using the standard minimum distance criterion of the central field model.
The analytical interacting hardsphere model cannot give the probability of sticking or rebounding after collisions. To assess the fraction of collisions leading to a stable product cluster, we analyzed the large data set of collisions obtained during the validation of the analytical interacting hardsphere model, taking into account the limitations of the classical force fields employed, as well as the timescales of equilibration processes in the atmosphere, not simulated explicitly in the binary collision setup.
Our analysis shows that if a monomer can stay on the cluster surface for more than a characteristic time (typically tens of picoseconds) after a collision, then a stable product cluster is formed. With this criterion, the sticking probability, sticking rate coefficient, and mass accommodation coefficient can be calculated by analyzing a set of molecular dynamics collision trajectories where impact parameters and relative speeds are properly sampled.
The mass accommodation coefficient decreases with the increasing temperature and decreasing cluster size. For collisions between H_{2}SO_{4} and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n}, the unit mass accommodation coefficient assumption is generally valid in an atmospherically relevant temperature range of 200–400 K, regardless of the cluster size. However, the assumption may not hold if the temperature is further increased. For collisions between (CH_{3})_{2}NH and [${\mathrm{HSO}}_{\mathrm{4}}^{}$ ⋅ $({\mathrm{CH}}_{\mathrm{3}}{)}_{\mathrm{2}}{\mathrm{NH}}_{\mathrm{2}}^{+}$]_{n}, the unit mass accommodation coefficient assumption is approximately valid at 200 K but becomes significantly worse at elevated temperatures; e.g., at 400 K it ranges from 0.1–0.6 for clusters with 1–32 monomers.
Future research directions include the investigation of systems where the transition between the two regimes defining the critical impact parameter occurs at lower and more relevant relative initial velocities, as well as the application of the analytical interacting hardsphere model to cluster–cluster collisions. The model presented in this paper is expected to provide a useful tool for the atmospheric NPF community due to its relative simplicity, demonstrated good accuracy, and widespread applicability.
The data generated in this study are fully represented in the figures and tables shown in the paper and the Supplement. Input files for simulations are available from the authors upon reasonable request.
The supplement related to this article is available online at: https://doi.org/10.5194/acp2359932023supplement.
HY conceived the research and derived the theoretical framework of the model. IN and HY performed the molecular dynamics collision simulations. IN provided theoretical background for the central field approach. BR and VT performed the metadynamics simulations. VT and JK helped to prepare the clusters used in the simulation. TK and HV helped to plan the project. HY, IN, and BR analyzed the simulation data and wrote the manuscript with comments from all other authors.
The contact author has declared that none of the authors has any competing interests.
Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Computational resources were provided by the CSC–IT Center for Science Ltd., Finland. The authors wish to thank the Finnish Computing Competence Infrastructure (FCCI) for supporting this project with computational and data storage resources.
This research has been supported by the European Research Council (project no. 692891
DAMOCLES), the Academy of Finland Flagship Program (grant no. 337549 ACCC) and
Centers of Excellence Program (grant no. 346368 VILMA).
Openaccess funding was provided by the Helsinki University Library.
This paper was edited by Hinrich Grothe and reviewed by two anonymous referees.
Barducci, A., Bussi, G., and Parrinello, M.: WellTempered Metadynamics: A Smoothly Converging and Tunable FreeEnergy Method, Phys. Rev. Lett., 100, 020603, https://doi.org/10.1103/PhysRevLett.100.020603, 2008. a
Becker, R. and Döring, W.: Kinetische Behandlung der Keimbildung in übersättigten Dämpfen, Ann. Phys., 416, 719–752, https://doi.org/10.1002/andp.19354160806, 1935. a
Clary, D.: Calculations of rate constants for ionmolecule reactions using a combined capture and centrifugal sudden approximation, Mol. Phys., 54, 605–618, https://doi.org/10.1080/00268978500100461, 1985. a
Elm, J., Jen, C. N., Kurtén, T., and Vehkamäki, H.: Strong Hydrogen Bonded Molecular Interactions between Atmospheric Diamines and Sulfuric Acid, J. Phys. Chem. A, 120, 3693–3700, https://doi.org/10.1021/acs.jpca.6b03192, 2016. a
Farkas, L.: Keimbildungsgeschwindigkeit in übersättigten Dämpfen, Z. Phys. Chem., 125U, 236–242, https://doi.org/10.1515/zpch192712513, 1927. a
Fuchs, N. and Sutugin, A.: Coagulation rate of highly dispersed aerosols, J. Colloid Sci., 20, 492–500, https://doi.org/10.1016/00958522(65)900310, 1965. a, b
Gioumousis, G. and Stevenson, D. P.: Reactions of Gaseous Molecule Ions with Gaseous Molecules. V. Theory, J. Chem. Phys., 29, 294–299, https://doi.org/10.1063/1.1744477, 1958. a
Guo, S., Hu, M., Zamora, M. L., Peng, J., Shang, D., Zheng, J., Du, Z., Wu, Z., Shao, M., Zeng, L., Molina, M. J., and Zhang, R.: Elucidating severe urban haze formation in China, P. Natl. Acad. Sci. USA, 111, 17373–17378, https://doi.org/10.1073/pnas.1419604111, 2014. a
Hallquist, M., Wenger, J. C., Baltensperger, U., Rudich, Y., Simpson, D., Claeys, M., Dommen, J., Donahue, N. M., George, C., Goldstein, A. H., Hamilton, J. F., Herrmann, H., Hoffmann, T., Iinuma, Y., Jang, M., Jenkin, M. E., Jimenez, J. L., KiendlerScharr, A., Maenhaut, W., McFiggans, G., Mentel, Th. F., Monod, A., Prévôt, A. S. H., Seinfeld, J. H., Surratt, J. D., Szmigielski, R., and Wildt, J.: The formation, properties and impact of secondary organic aerosol: current and emerging issues, Atmos. Chem. Phys., 9, 5155–5236, https://doi.org/10.5194/acp951552009, 2009. a
Halonen, R., Zapadinsky, E., Kurtén, T., Vehkamäki, H., and Reischl, B.: Rate enhancement in collisions of sulfuric acid molecules due to longrange intermolecular forces, Atmos. Chem. Phys., 19, 13355–13366, https://doi.org/10.5194/acp19133552019, 2019. a, b, c, d
Hamaker, H.: The London—van der Waals attraction between spherical particles, Physica, 4, 1058–1072, https://doi.org/10.1016/S00318914(37)802037, 1937. a, b
Jorgensen, W. L., Maxwell, D. S., and TiradoRives, J.: Development and Testing of the OPLS AllAtom Force Field on Conformational Energetics and Properties of Organic Liquids, J. Am. Chem. Soc., 118, 11225–11236, https://doi.org/10.1021/ja9621760, 1996. a
Kirkby, J., Curtius, J., Almeida, J., Dunne, E., Duplissy, J., Ehrhart, S., Franchin, A., Gagné, S., Ickes, L., Kürten, A., Kupc, A., Metzger, A., Riccobono, F., Rondo, L., Schobesberger, S., Tsagkogeorgas, G., Wimmer, D., Amorim, A., Bianchi, F., Breitenlechner, M., David, A., Dommen, J., Downard, A., Ehn, M., Flagan, R. C., Haider, S., Hansel, A., Hauser, D., Jud, W., Junninen, H., Kreissl, F., Kvashin, A., Laaksonen, A., Lehtipalo, K., Lima, J., Lovejoy, E. R., Makhmutov, V., Mathot, S., Mikkilä, J., Minginette, P., Mogo, S., Nieminen, T., Onnela, A., Pereira, P., Petäjä, T., Schnitzhofer, R., Seinfeld, J. H., Mikko Sipilä, Y. S., Stratmann, F., Tomé, A., Vanhanen, J., Viisanen, Y., Vrtala, A., Wagner, P. E., Walther, H., Weingartner, E., Wex, H., Winkler, P. M., Carslaw, K. S., Worsnop, D. R., Baltensperger, U., and Kulmala, M.: Role of sulphuric acid, ammonia and galactic cosmic rays in atmospheric aerosol nucleation, Nature, 476, 429–433, https://doi.org/10.1038/nature10343, 2011. a
Kulmala, M., Pirjola, L., and Mäkelä, J. M.: Elucidating severe urban haze formation in China, P. Natl. Acad. Sci. USA, 404, 66–69, https://doi.org/10.1038/35003550, 2000. a
Kulmala, M., Kontkanen, J., Junninen, H., Lehtipalo, K., Manninen, H. E., Nieminen, T., Petäjä, T., Sipilä, M., Schobesberger, S., Rantala, P., Franchin, A., Jokinen, T., Järvinen, E., Äijälä, M., Kangasluoma, J., Hakala, J., Aalto, P. P., Paasonen, P., Mikkilä, J., Vanhanen, J., Aalto, J., Hakola, H., Makkonen, U., Ruuskanen, T., Mauldin, R. L., Duplissy, J., Vehkamäki, H., Bäck, J., Kortelainen, A., Riipinen, I., Kurtén, T., Johnston, M. V., Smith, J. N., Ehn, M., Mentel, T. F., Lehtinen, K. E. J., Laaksonen, A., Kerminen, V.M., and Worsnop, D. R.: Direct Observations of Atmospheric Aerosol Nucleation, Science, 339, 943–946, https://doi.org/10.1126/science.1227385, 2013. a
Kurtén, T., Loukonen, V., Vehkamäki, H., and Kulmala, M.: Amines are likely to enhance neutral and ioninduced sulfuric acidwater nucleation in the atmosphere more effectively than ammonia, Atmos. Chem. Phys., 8, 4095–4103, https://doi.org/10.5194/acp840952008, 2008. a
Landau, L. D. and Lifshitz, E. M.: Mechanics, vol. 1, ButterworthHeinemann, https://doi.org/10.1016/C20090255693, 1976. a
Lehtipalo, K., Rondo, L., Kontkanen, J., Schobesberger, S., Jokinen, T., Sarnela, N., Kürten, A., Ehrhart, S., Franchin, A., Nieminen, T., Riccobono, F., Sipilä, M., YliJuuti, T., Duplissy, J., Adamov, A., Ahlm, L., Almeida, J., Amorim, A., Bianchi, F., Breitenlechner, M., Dommen, J., Downard, A. J., Dunne, E. M., Flagan, R. C., Guida, R., Hakala, J., Hansel, A., Jud, W., Kangasluoma, J., Kerminen, V.M., Keskinen, H., Kim, J., Kirkby, J., Kupc, A., KupiainenMäättä, O., Laaksonen, A., Lawler, M. J., Leiminger, M., Mathot, S., Olenius, T., Ortega, I. K., Onnela, A., Petäjä, T., Praplan, A., Rissanen, M. P., Ruuskanen, T., Santos, F. D., Schallhart, S., Schnitzhofer, R., Simon, M., Smith, J. N., Tröstl, J., Tsagkogeorgas, G., Tomé, A., Vaattovaara, P., Vehkamäki, H., Vrtala, A. E., Wagner, P. E., Williamson, C., Wimmer, D., Winkler, P. M., Virtanen, A., Donahue, N. M., Carslaw, K. S., Baltensperger, U., Riipinen, I., Curtius, J., Worsnop, D. R., and Kulmala, M.: The effect of acid–base clustering and ions on the growth of atmospheric nanoparticles, Nat. Commun., 7, 11594, https://doi.org/10.1038/ncomms11594, 2016. a
Leite, F. L., Bueno, C. C., Da Róz, A. L., Ziemath, E. C., and Oliveira, O. N.: Theoretical Models for Surface Forces and Adhesion and Their Measurement Using Atomic Force Microscopy, Int. J. Mol. Sci., 13, 12773–12856, https://doi.org/10.3390/ijms131012773, 2012. a
Loukonen, V., Kurtén, T., Ortega, I. K., Vehkamäki, H., Pádua, A. A. H., Sellegri, K., and Kulmala, M.: Enhancing effect of dimethylamine in sulfuric acid nucleation in the presence of water – a computational study, Atmos. Chem. Phys., 10, 4961–4974, https://doi.org/10.5194/acp1049612010, 2010. a, b
Moran, T. F. and Hamill, W. H.: Cross Sections of Ion–PermanentDipole Reactions by Mass Spectrometry, J. Chem. Phys., 39, 1413–1422, https://doi.org/10.1063/1.1734457, 1963. a
Myllys, N., Olenius, T., Kurtén, T., Vehkamäki, H., Riipinen, I., and Elm, J.: Effect of Bisulfate, Ammonia, and Ammonium on the Clustering of Organic Acids and Sulfuric Acid, J. Phys. Chem. A, 121, 4812–4824, https://doi.org/10.1021/acs.jpca.7b03981, 2017. a
Neefjes, I., Halonen, R., Vehkamäki, H., and Reischl, B.: Modeling approaches for atmospheric ion–dipole collisions: allatom trajectory simulations and central field methods, Atmos. Chem. Phys., 22, 11155–11172, https://doi.org/10.5194/acp22111552022, 2022. a, b, c
Ortega, I. K., Kupiainen, O., Kurtén, T., Olenius, T., Wilkman, O., McGrath, M. J., Loukonen, V., and Vehkamäki, H.: From quantum chemical formation free energies to evaporation rates, Atmos. Chem. Phys., 12, 225–235, https://doi.org/10.5194/acp122252012, 2012. a
Ouyang, H., Gopalakrishnan, R., and Hogan, C. J.: Nanoparticle collisions in the gas phase in the presence of singular contact potentials, J. Chem. Phys., 137, 064316, https://doi.org/10.1063/1.4742064, 2012. a
Plimpton, S.: Fast Parallel Algorithms for ShortRange Molecular Dynamics, J. Comput. Phys., 117, 1–19, https://doi.org/10.1006/jcph.1995.1039, 1995. a, b
Pope III, C. A. and Dockery, D. W.: Aerosols, Climate, and the Hydrological Cycle, J. Air Waste Manage., 56, 709–742, https://doi.org/10.1080/10473289.2006.10464485, 2006. a
Pöschl, U.: Atmospheric Aerosols: Composition, Transformation, Climate and Health Effects, Angew. Chem. Int. Edit., 44, 7520–7540, https://doi.org/10.1002/anie.200501122, 2005. a
Ramanathan, V., Crutzen, P. J., Kiehl, J. T., and Rosenfeld, D.: Aerosols, Climate, and the Hydrological Cycle, Science, 294, 2119–2124, https://doi.org/10.1126/science.1064034, 2001. a
Sceats, M. G.: Brownian coagulation in aerosols—the role of long range forces, J. Colloid Interface Sci., 129, 105–112, https://doi.org/10.1016/00219797(89)904190, 1989. a
Schenter, G. K., Kathmann, S. M., and Garrett, B. C.: Dynamical Nucleation Theory: A New Molecular Approach to VaporLiquid Nucleation, Phys. Rev. Lett., 82, 3484–3487, https://doi.org/10.1103/PhysRevLett.82.3484, 1999. a
Seinfeld, J. H., Bretherton, C., Carslaw, K. S., Coe, H., DeMott, P. J., Dunlea, E. J., Feingold, G., Ghan, S., Guenther, A. B., Kahn, R., Kraucunas, I., Kreidenweis, S. M., Molina, M. J., Nenes, A., Penner, J. E., Prather, K. A., Ramanathan, V., Ramaswamy, V., Rasch, P. J., Ravishankara, A. R., Rosenfeld, D., Stephens, G., and Wood, R.: Improving our fundamental understanding of the role of aerosolcloud interactions in the climate system, P. Natl. Acad. Sci. USA, 113, 5781–5790, https://doi.org/10.1073/pnas.1514043113, 2016. a
Sipilä, M., Berndt, T., Petäjä, T., Brus, D., Vanhanen, J., Stratmann, F., Patokoski, J., Mauldin, R. L., Hyvärinen, A.P., Lihavainen, H., and Kulmala, M.: The Role of Sulfuric Acid in Atmospheric Nucleation, Science, 327, 1243–1246, https://doi.org/10.1126/science.1180315, 2010. a
Stolzenburg, D., Simon, M., Ranjithkumar, A., Kürten, A., Lehtipalo, K., Gordon, H., Ehrhart, S., Finkenzeller, H., Pichelstorfer, L., Nieminen, T., He, X.C., Brilke, S., Xiao, M., Amorim, A., Baalbaki, R., Baccarini, A., Beck, L., Bräkling, S., Caudillo Murillo, L., Chen, D., Chu, B., Dada, L., Dias, A., Dommen, J., Duplissy, J., El Haddad, I., Fischer, L., Gonzalez Carracedo, L., Heinritzi, M., Kim, C., Koenig, T. K., Kong, W., Lamkaddam, H., Lee, C. P., Leiminger, M., Li, Z., Makhmutov, V., Manninen, H. E., Marie, G., Marten, R., Müller, T., Nie, W., Partoll, E., Petäjä, T., Pfeifer, J., Philippov, M., Rissanen, M. P., Rörup, B., Schobesberger, S., Schuchmann, S., Shen, J., Sipilä, M., Steiner, G., Stozhkov, Y., Tauber, C., Tham, Y. J., Tomé, A., VazquezPufleau, M., Wagner, A. C., Wang, M., Wang, Y., Weber, S. K., Wimmer, D., Wlasits, P. J., Wu, Y., Ye, Q., ZaunerWieczorek, M., Baltensperger, U., Carslaw, K. S., Curtius, J., Donahue, N. M., Flagan, R. C., Hansel, A., Kulmala, M., Lelieveld, J., Volkamer, R., Kirkby, J., and Winkler, P. M.: Enhanced growth rate of atmospheric particles from sulfuric acid, Atmos. Chem. Phys., 20, 7359–7372, https://doi.org/10.5194/acp2073592020, 2020. a
Su, T. and Bowers, M. T.: Theory of ion‐polar molecule collisions. Comparison with experimental charge transfer reactions of rare gas ions to geometric isomers of difluorobenzene and dichloroethylene, J. Chem. Phys., 58, 3027–3037, https://doi.org/10.1063/1.1679615, 1973. a
Su, T., Su, E. C., and Bowers, M. T.: Ion–polar molecule collisions. Conservation of angular momentum in the average dipole orientation theory. The AADO theory, J. Chem. Phys., 69, 2243–2250, https://doi.org/10.1063/1.436783, 1978. a
Temelso, B., Phan, T. N., and Shields, G. C.: Computational Study of the Hydration of Sulfuric Acid Dimers: Implications for Acid Dissociation and Aerosol Formation, J. Phys. Chem. A, 116, 9745–9758, https://doi.org/10.1021/jp3054394, 2012. a
Tribello, G. A., Bonomi, M., Branduardi, D., Camilloni, C., and Bussi, G.: PLUMED 2: New feathers for an old bird, Comput. Phys. Commun., 185, 604–613, https://doi.org/10.1016/j.cpc.2013.09.018, 2014. a
Vehkamäki, H. and Riipinen, I.: Thermodynamics and kinetics of atmospheric aerosol particle formation and growth, Chem. Soc. Rev., 41, 5160–5173, https://doi.org/10.1039/C2CS00002D, 2012. a
Yang, H., Goudeli, E., and Hogan, C. J.: Condensation and dissociation rates for gas phase metal clusters from molecular dynamics trajectory calculations, J. Chem. Phys., 148, 164304, https://doi.org/10.1063/1.5026689, 2018. a, b
Yang, H., Drossinos, Y., and Hogan, C. J.: Excess thermal energy and latent heat in nanocluster collisional growth, J. Chem. Phys., 151, 224304, https://doi.org/10.1063/1.5129918, 2019. a
Yang, H., Song, G., and Hogan, C. J.: A molecular dynamics study of collisional heat transfer to nanoclusters in the gas phase, J. Aerosol Sci., 159, 105891, https://doi.org/10.1016/j.jaerosci.2021.105891, 2022. a, b
Zhang, R., Khalizov, A., Wang, L., Hu, M., and Xu, W.: Nucleation and Growth of Nanoparticles in the Atmosphere, Chem. Rev., 112, 1957–2011, https://doi.org/10.1021/cr2001756, 2012. a