the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Gas-phase collision rate enhancement factors for acid-base clusters up to 2 nm in diameter from atomistic simulation and the interacting hard sphere model
Abstract. Collisions of neutral molecules and clusters is the most prevalent pathway in atmospheric new particle formation (NPF), with direct implications on air quality and climate. Until recently, these collisions have been modeled mainly using non-interacting hard sphere (NHS) models, which systematically underestimate collision and particle formation rates due to omission of long-range interactions. Lately, atomistic simulations which account for long-range interactions have been used to study neutral molecule-molecule and molecule-cluster collisions, but studies on cluster-cluster collisions have still been lacking despite the relevant role they can play e.g. in haze formation in polluted urban areas. We have therefore studied collisions between neutral clusters of N bisulphate and N dimethylammonium ions at T = 300 K up to N = 32 using atomistic molecular dynamics (MD) simulations. Direct simulation results have then been compared against both the traditional NHS model and the newly proposed interacting hard sphere (IHS) variant. We find the collision rates given by the NHS to be enhanced by a factor of 2.18–5.61 in the atomistic MD simulations, with enhancement decreasing with cluster size, and an asymptotic limit ≈ 2. The IHS model yields a constant enhancement factor of 3.36 over the NHS model for collisions between same-sized clusters, which decreases with increasing cluster size ratio. Our results demonstrate how even collisions between clusters of tens of acid-base pairs at a relatively high temperature cannot be accurately modeled by neglecting long-range interactions. We also show that the MD results cannot be reproduced by simple point-particle models, highlighting the importance of atomistic details of intermolecular interactions.
- Preprint
(4446 KB) - Metadata XML
- BibTeX
- EndNote
Status: closed
-
RC1: 'Comment on egusphere-2025-507', Anonymous Referee #1, 24 Mar 2025
Tikkanen et al. studied sulfuric acid-dimethyalmine cluster-cluster collisions. The hard-sphere (HS) model was compared to the interactive hard-sphere (IHS) model using the Hamaker approach, and benchmarked against direct trajectory simulations using the OPLS force-field. The authors find that both the direct and IHS model predict larger collision coefficients as they do not neglect long-range interactions, suggesting that the commonly used models should predict larger new particle formation rates. The authors also find that the direct and IHS model do not agree, due to problems with fitting the potential of mean force, and the underlying assumptions in the Hamaker approach, suggesting that for complex systems, the direct approach is needed.
The paper is well written, and the science is sound; I support its publication in Atmospheric Chemistry and Physics, but I have a few minor comments that will help clarify the authors’ intent.
Minor comments:
Throughout the paper [HSO4- (CH3)2NH2+] is sometimes referred to as a monomer. I understand that this is because it is treated as the “monomer” in the IHS Hamaker approach, but this terminology can easily lead to confusion. In the field of atmospheric science, this species is typically defined as a dimer. Calling it a “dimer unit” would clarify the meaning without requiring extensive changes to the surrounding text. Alternatively, the authors could explicitly state this naming convention early in the introduction.
Page 6 – line 171: For larger SA-DMA clusters, all recent quantum mechanical calculations indicate that clusters with complete proton transfer yield the lowest energy. So I do not think it is a ”questionable” assumption at all.
Page 8 – line 209: Can the authors explain why the constraint on the radius of gyration was applied? Is the dimer unit unstable with the specific force-field used?
Page 9 – Section 3.2: The section focuses heavily on the “attractive tail”; however, the tail region is not explicitly defined. From the plots, I assume it is defined from the shoulder at the largest r value and beyond.
Technical comment:
Page 6 – line 146: clusterd → cluster
Page 11 – line 265: Perhaps the authors could specify that the different numbers of MD trajectories were due to computational constraints.
Page 12 – Figure 4: Could the 2RHs values be added to the plot to make it easier to follow?
Figure 5, 7, 8, and 9: Add that N refers to number of the SA-DMA units: [HSO4- (CH3)2NH2+]N , to ensure that the figure are self-explanatory.
Citation: https://doi.org/10.5194/egusphere-2025-507-RC1 -
RC2: 'Comment on egusphere-2025-507', Anonymous Referee #2, 03 Apr 2025
This study investigates the collision dynamics of neutral bisulfate-dimethylammonium clusters using atomistic MD simulations. By comparing with Non-interacting Hard Sphere (NHS) and Interacting Hard Sphere (IHS) models, it quantifies the enhancement of collision rates due to long-range interactions (e.g., dipole-dipole, van der Waals). Current simulations like ACDC models rely on empirical enhancement factors (e.g., βij(N)), but lack data for large clusters (N > 20). This work provides a benchmark for refining size-dependent parameterizations, at the same time, The IHS model’s constant enhancement factor (WIHS=3.36) fails to capture the size-dependent decay observed in MD (asymptotic limit ≈ 2 for N=32), highlighting its limitations for large clusters. Thus, It demonstrates that simplified point-particle models (like IHS) cannot fully replace atomistic simulations, especially in polluted regions where cluster complexity (e.g., surface mobility, many-body effects) dominates. This manuscript is clearly articulated and presents robust scientific findings; I endorse its publication in Atmospheric Chemistry and Physics, before the following issues were further clarified:
- Page 3, line 64-65, I note that only collision between clusters of the same size are studied. Atmospheric NPF involves clusters of diverse sizes (N=1 to N>100N), the same-size collisions (Ni=Nj) would ignore dominant asymmetric cases, missing critical dynamics (e.g., dipole-induced dipole forces in R≫1collisions). potentially underestimating rates. The limitations should be discussed here.
- Page 5, line 139-140, The processes (numerical techniques) to find the minimum of the function ω(r) seems to not clearly Additional description need to be provided, asthe function ω(r) is the key parameter to calculate the enhancement factor of the IHS model.
- Page 7, line 195, the authors mentioned that the value of Rb = 0.2RHS was ultimately used, after they test different buffer values Rb, and found the results were not sensitive to this choice. Still, it not clear for the using the value of Rb = 0.2RHS.
Citation: https://doi.org/10.5194/egusphere-2025-507-RC2 -
AC1: 'Comment on egusphere-2025-507', Valtteri Tikkanen, 03 Jun 2025
We thank both referees for their positive and constructive feedback, which we carefully considered in the revised manuscript. In our opinion, all the referee comments have been taken into account in a satisfactory manner, as listed in detail below.
RC1
Minor comments
- Throughout the paper [HSO4-(CH3)2NH2+]is sometimes referred to as a monomer. I understand that this is because it is treated as the “monomer” in the IHS Hamaker approach, but this terminology can easily lead to confusion. In the field of atmospheric science, this species is typically defined as a dimer. Calling it a “dimer unit” would clarify the meaning without requiring extensive changes to the surrounding text. Alternatively, the authors could explicitly state this naming convention early in the introduction.
We agree that the concept of a ”monomer” might be misunderstood by some readers. Therefore, in the revised text the meaning of a monomer is given explicitly: ” In this work, the term `monomer' refers to the acid-base dimer unit, [HSO4- (CH3)2NH2+]1 , because this `heterodimer' is the logical unit in modeling and simulating neutral bisulphate dimethylammonium clusters.”
- For larger SA-DMA clusters, all recent quantum mechanical calculations indicate that clusters with complete proton transfer yield the lowest energy. So I do not think it is a ”questionable” assumption at all.
The comment about full proton transfer being ”questionable” for these systems has been removed, as suggested, and the sentence modfied:
“For clusters N > 4, we have also assumed complete proton transfer, i.e., clusters consisting only of bisulphate and dimethylammonium ions, in line with the QM minimum energy structures of the largest clusters available in the database (Elm, 2019).”
- Can the authors explain why the constraint on the radius of gyration was applied? Is the dimer unit unstable with the specific force-field used?
The dimer unit of bisulphate + dimethylammonium is indeed very stable in the OPLS-AA forcefield description. We have clarified the reasoning behind constraining the radius of gyration as follows: ” The radii of gyration of each [HSO4- (CH3)2NH2+]N monomer were constrained by a harmonic upper wall with a spring constant of k = 10 eV/Å2, starting at the value of Rg = 2.5 Å, to ensure the dimer units remain intact at intermediate distances, while still allowing for necessary rearrangements upon cluster formation.”
Also, we added the following text to the discussion regarding larger clusters (for which this constrain becomes more important): ” For these simulations, harmonic upper walls were put on the value of the radii of gyration at distances deemed appropriate to achieve a compromise between avoiding elongation of the clusters leading to coalescence already at intermediate distances, and too rigidly constraining the clusters to their original spherical geometries.”
- Section 3.2: The section focuses heavily on the “attractive tail”; however, the tail region is not explicitly defined. From the plots, I assume it is defined from the shoulder at the largestrvalue and beyond.
This is correct, the meaning of the tail region is now explicitly defined:
“For N = 1, the PMF exhibits a global minimum at r = 4.05 Å with a well depth of 1.53 eV, and the attractive tail for r > 8 Å is in excellent agreement with the rotationally-averaged dipole-dipole interaction (Eq. 16) using the average dipole moment from simulation (see Table 1), indicating that the long-range attractive interaction between the heterodimers can be effectively modeled as that between two point dipoles.”
and
“It must be noted that for systems N > 1 we solely focus on the tails of the PMF curves, i.e. distances at which the PMF increases monotonically to zero, for multiple reasons: […]”
Technical comments
- Page 6 – line 146: clusterd → cluster: has been corrected
- Page 11 – line 265: Perhaps the authors could specify that the different numbers of MD trajectories were due to computational constraints.
We have modified the sentence:
“We performed ntraj = 500, 500, 200, 200 and 100 individual MD trajectory simulations for collisions between N = 1,2,4,8 and 32 [HSO4- (CH3)2NH2+]N , clusters, respectively, to balance between accuracy and increasing computational cost for larger system sizes.”
- Figure 4: Could the 2RHsvalues be added to the plot to make it easier to follow?
A dashed grey line has been added at r=2RHS to all panels in the revised Figure 4.
- Figure 5, 7, 8, and 9: Add thatNrefers to number of the SA-DMA units: [HSO4- (CH3)2NH2+]N , to ensure that the figure are self-explanatory.
We have modified the captions of Figures 5, 7, 8 and 9 accordingly.
RC2
- Page 3, line 64-65, I note that only collision between clusters of the same size are studied. Atmospheric NPF involves clusters of diverse sizes (N=1 to N>100N), the same-size collisions (Ni=Nj) would ignore dominant asymmetric cases, missing critical dynamics (e.g., dipole-induced dipole forces in R≫1collisions). potentially underestimating rates. The limitations should be discussed here.
We agree that for the full picture of NPF in the atmosphere, collisions between very diverse molecules and clusters need to be considered. Here, we focus on the same-sized cluster-cluster systems as a continuation to our previous work considering molecule-molecule and molecule-cluster collisions. To the revised text we have added the following clarification:
”While the limitation to same-sized cluster collisions in the MD simulations in this work does not account for the vast majority of asymmetric collisions in real atmospheric processes, it still provides a useful starting point to investigate size-dependent collision rate enhancements, when attractive cluster-cluster interactions are taken into account.”
- Page 5, line 139-140, The processes (numerical techniques) to find the minimum of the function ω(r) seems to not clearly Additional description need to be provided, as the function ω(r) is the key parameter to calculate the enhancement factor of the IHS model.
The revised text has been expanded to make the procedure easier to follow:
”Consequently, we first determine the minimum of ω(r) = r2(1 − 2Ucc(r)/µv20) numerically. This numerical minimum is subsequently substituted into Eq. 5, facilitating the determination of to determine the critical impact parameter. The obtained critical impact parameter,
b2c= ωv (rmin), if rmin > 2RHS,
or b2c=ωv (RHS), if rmin ≤ 2RHS, (11)
is then used in Eqs. 6 and 7 to calculate the enhancement factor of the IHS model numerically.”
- Page 7, line 195, the authors mentioned that the value of Rb = 0.2RHS was ultimately used, after they test different buffer values Rb, and found the results were not sensitive to this choice. Still, it not clear for the using the value of Rb = 0.2RHS.
The buffer value Rb was first introduced to avoid missing successful collisions in the analysis of the MD collision trajectories, in case clusters changed shape before collision, which could potentially lead to center of mass distances exceeding the sum of hard-sphere radii in the product cluster. However, the results are not sensitive to this choice, and for typical collisions, the center of mass distances were actually significantly smaller than 2RHS, as can be seen in Figure 8 in section 3.3, which makes this point clearer.
We have added the following sentence to the revised manuscript:
“In fact, for typical collisions r < 2RHS (see Fig. 8a in Sec. 3.3).”
Citation: https://doi.org/10.5194/egusphere-2025-507-AC1
Status: closed
-
RC1: 'Comment on egusphere-2025-507', Anonymous Referee #1, 24 Mar 2025
Tikkanen et al. studied sulfuric acid-dimethyalmine cluster-cluster collisions. The hard-sphere (HS) model was compared to the interactive hard-sphere (IHS) model using the Hamaker approach, and benchmarked against direct trajectory simulations using the OPLS force-field. The authors find that both the direct and IHS model predict larger collision coefficients as they do not neglect long-range interactions, suggesting that the commonly used models should predict larger new particle formation rates. The authors also find that the direct and IHS model do not agree, due to problems with fitting the potential of mean force, and the underlying assumptions in the Hamaker approach, suggesting that for complex systems, the direct approach is needed.
The paper is well written, and the science is sound; I support its publication in Atmospheric Chemistry and Physics, but I have a few minor comments that will help clarify the authors’ intent.
Minor comments:
Throughout the paper [HSO4- (CH3)2NH2+] is sometimes referred to as a monomer. I understand that this is because it is treated as the “monomer” in the IHS Hamaker approach, but this terminology can easily lead to confusion. In the field of atmospheric science, this species is typically defined as a dimer. Calling it a “dimer unit” would clarify the meaning without requiring extensive changes to the surrounding text. Alternatively, the authors could explicitly state this naming convention early in the introduction.
Page 6 – line 171: For larger SA-DMA clusters, all recent quantum mechanical calculations indicate that clusters with complete proton transfer yield the lowest energy. So I do not think it is a ”questionable” assumption at all.
Page 8 – line 209: Can the authors explain why the constraint on the radius of gyration was applied? Is the dimer unit unstable with the specific force-field used?
Page 9 – Section 3.2: The section focuses heavily on the “attractive tail”; however, the tail region is not explicitly defined. From the plots, I assume it is defined from the shoulder at the largest r value and beyond.
Technical comment:
Page 6 – line 146: clusterd → cluster
Page 11 – line 265: Perhaps the authors could specify that the different numbers of MD trajectories were due to computational constraints.
Page 12 – Figure 4: Could the 2RHs values be added to the plot to make it easier to follow?
Figure 5, 7, 8, and 9: Add that N refers to number of the SA-DMA units: [HSO4- (CH3)2NH2+]N , to ensure that the figure are self-explanatory.
Citation: https://doi.org/10.5194/egusphere-2025-507-RC1 -
RC2: 'Comment on egusphere-2025-507', Anonymous Referee #2, 03 Apr 2025
This study investigates the collision dynamics of neutral bisulfate-dimethylammonium clusters using atomistic MD simulations. By comparing with Non-interacting Hard Sphere (NHS) and Interacting Hard Sphere (IHS) models, it quantifies the enhancement of collision rates due to long-range interactions (e.g., dipole-dipole, van der Waals). Current simulations like ACDC models rely on empirical enhancement factors (e.g., βij(N)), but lack data for large clusters (N > 20). This work provides a benchmark for refining size-dependent parameterizations, at the same time, The IHS model’s constant enhancement factor (WIHS=3.36) fails to capture the size-dependent decay observed in MD (asymptotic limit ≈ 2 for N=32), highlighting its limitations for large clusters. Thus, It demonstrates that simplified point-particle models (like IHS) cannot fully replace atomistic simulations, especially in polluted regions where cluster complexity (e.g., surface mobility, many-body effects) dominates. This manuscript is clearly articulated and presents robust scientific findings; I endorse its publication in Atmospheric Chemistry and Physics, before the following issues were further clarified:
- Page 3, line 64-65, I note that only collision between clusters of the same size are studied. Atmospheric NPF involves clusters of diverse sizes (N=1 to N>100N), the same-size collisions (Ni=Nj) would ignore dominant asymmetric cases, missing critical dynamics (e.g., dipole-induced dipole forces in R≫1collisions). potentially underestimating rates. The limitations should be discussed here.
- Page 5, line 139-140, The processes (numerical techniques) to find the minimum of the function ω(r) seems to not clearly Additional description need to be provided, asthe function ω(r) is the key parameter to calculate the enhancement factor of the IHS model.
- Page 7, line 195, the authors mentioned that the value of Rb = 0.2RHS was ultimately used, after they test different buffer values Rb, and found the results were not sensitive to this choice. Still, it not clear for the using the value of Rb = 0.2RHS.
Citation: https://doi.org/10.5194/egusphere-2025-507-RC2 -
AC1: 'Comment on egusphere-2025-507', Valtteri Tikkanen, 03 Jun 2025
We thank both referees for their positive and constructive feedback, which we carefully considered in the revised manuscript. In our opinion, all the referee comments have been taken into account in a satisfactory manner, as listed in detail below.
RC1
Minor comments
- Throughout the paper [HSO4-(CH3)2NH2+]is sometimes referred to as a monomer. I understand that this is because it is treated as the “monomer” in the IHS Hamaker approach, but this terminology can easily lead to confusion. In the field of atmospheric science, this species is typically defined as a dimer. Calling it a “dimer unit” would clarify the meaning without requiring extensive changes to the surrounding text. Alternatively, the authors could explicitly state this naming convention early in the introduction.
We agree that the concept of a ”monomer” might be misunderstood by some readers. Therefore, in the revised text the meaning of a monomer is given explicitly: ” In this work, the term `monomer' refers to the acid-base dimer unit, [HSO4- (CH3)2NH2+]1 , because this `heterodimer' is the logical unit in modeling and simulating neutral bisulphate dimethylammonium clusters.”
- For larger SA-DMA clusters, all recent quantum mechanical calculations indicate that clusters with complete proton transfer yield the lowest energy. So I do not think it is a ”questionable” assumption at all.
The comment about full proton transfer being ”questionable” for these systems has been removed, as suggested, and the sentence modfied:
“For clusters N > 4, we have also assumed complete proton transfer, i.e., clusters consisting only of bisulphate and dimethylammonium ions, in line with the QM minimum energy structures of the largest clusters available in the database (Elm, 2019).”
- Can the authors explain why the constraint on the radius of gyration was applied? Is the dimer unit unstable with the specific force-field used?
The dimer unit of bisulphate + dimethylammonium is indeed very stable in the OPLS-AA forcefield description. We have clarified the reasoning behind constraining the radius of gyration as follows: ” The radii of gyration of each [HSO4- (CH3)2NH2+]N monomer were constrained by a harmonic upper wall with a spring constant of k = 10 eV/Å2, starting at the value of Rg = 2.5 Å, to ensure the dimer units remain intact at intermediate distances, while still allowing for necessary rearrangements upon cluster formation.”
Also, we added the following text to the discussion regarding larger clusters (for which this constrain becomes more important): ” For these simulations, harmonic upper walls were put on the value of the radii of gyration at distances deemed appropriate to achieve a compromise between avoiding elongation of the clusters leading to coalescence already at intermediate distances, and too rigidly constraining the clusters to their original spherical geometries.”
- Section 3.2: The section focuses heavily on the “attractive tail”; however, the tail region is not explicitly defined. From the plots, I assume it is defined from the shoulder at the largestrvalue and beyond.
This is correct, the meaning of the tail region is now explicitly defined:
“For N = 1, the PMF exhibits a global minimum at r = 4.05 Å with a well depth of 1.53 eV, and the attractive tail for r > 8 Å is in excellent agreement with the rotationally-averaged dipole-dipole interaction (Eq. 16) using the average dipole moment from simulation (see Table 1), indicating that the long-range attractive interaction between the heterodimers can be effectively modeled as that between two point dipoles.”
and
“It must be noted that for systems N > 1 we solely focus on the tails of the PMF curves, i.e. distances at which the PMF increases monotonically to zero, for multiple reasons: […]”
Technical comments
- Page 6 – line 146: clusterd → cluster: has been corrected
- Page 11 – line 265: Perhaps the authors could specify that the different numbers of MD trajectories were due to computational constraints.
We have modified the sentence:
“We performed ntraj = 500, 500, 200, 200 and 100 individual MD trajectory simulations for collisions between N = 1,2,4,8 and 32 [HSO4- (CH3)2NH2+]N , clusters, respectively, to balance between accuracy and increasing computational cost for larger system sizes.”
- Figure 4: Could the 2RHsvalues be added to the plot to make it easier to follow?
A dashed grey line has been added at r=2RHS to all panels in the revised Figure 4.
- Figure 5, 7, 8, and 9: Add thatNrefers to number of the SA-DMA units: [HSO4- (CH3)2NH2+]N , to ensure that the figure are self-explanatory.
We have modified the captions of Figures 5, 7, 8 and 9 accordingly.
RC2
- Page 3, line 64-65, I note that only collision between clusters of the same size are studied. Atmospheric NPF involves clusters of diverse sizes (N=1 to N>100N), the same-size collisions (Ni=Nj) would ignore dominant asymmetric cases, missing critical dynamics (e.g., dipole-induced dipole forces in R≫1collisions). potentially underestimating rates. The limitations should be discussed here.
We agree that for the full picture of NPF in the atmosphere, collisions between very diverse molecules and clusters need to be considered. Here, we focus on the same-sized cluster-cluster systems as a continuation to our previous work considering molecule-molecule and molecule-cluster collisions. To the revised text we have added the following clarification:
”While the limitation to same-sized cluster collisions in the MD simulations in this work does not account for the vast majority of asymmetric collisions in real atmospheric processes, it still provides a useful starting point to investigate size-dependent collision rate enhancements, when attractive cluster-cluster interactions are taken into account.”
- Page 5, line 139-140, The processes (numerical techniques) to find the minimum of the function ω(r) seems to not clearly Additional description need to be provided, as the function ω(r) is the key parameter to calculate the enhancement factor of the IHS model.
The revised text has been expanded to make the procedure easier to follow:
”Consequently, we first determine the minimum of ω(r) = r2(1 − 2Ucc(r)/µv20) numerically. This numerical minimum is subsequently substituted into Eq. 5, facilitating the determination of to determine the critical impact parameter. The obtained critical impact parameter,
b2c= ωv (rmin), if rmin > 2RHS,
or b2c=ωv (RHS), if rmin ≤ 2RHS, (11)
is then used in Eqs. 6 and 7 to calculate the enhancement factor of the IHS model numerically.”
- Page 7, line 195, the authors mentioned that the value of Rb = 0.2RHS was ultimately used, after they test different buffer values Rb, and found the results were not sensitive to this choice. Still, it not clear for the using the value of Rb = 0.2RHS.
The buffer value Rb was first introduced to avoid missing successful collisions in the analysis of the MD collision trajectories, in case clusters changed shape before collision, which could potentially lead to center of mass distances exceeding the sum of hard-sphere radii in the product cluster. However, the results are not sensitive to this choice, and for typical collisions, the center of mass distances were actually significantly smaller than 2RHS, as can be seen in Figure 8 in section 3.3, which makes this point clearer.
We have added the following sentence to the revised manuscript:
“In fact, for typical collisions r < 2RHS (see Fig. 8a in Sec. 3.3).”
Citation: https://doi.org/10.5194/egusphere-2025-507-AC1
Viewed
HTML | XML | Total | BibTeX | EndNote | |
---|---|---|---|---|---|
399 | 133 | 18 | 550 | 20 | 40 |
- HTML: 399
- PDF: 133
- XML: 18
- Total: 550
- BibTeX: 20
- EndNote: 40
Viewed (geographical distribution)
Country | # | Views | % |
---|
Total: | 0 |
HTML: | 0 |
PDF: | 0 |
XML: | 0 |
- 1