Eva
Zunzunegui-Bru‡
a,
Elisabeth
Gruber
b,
Stefan
Bergmeister
b,
Miriam
Meyer
b,
Fabio
Zappa
b,
Massimiliano
Bartolomei
a,
Fernando
Pirani
c,
Pablo
Villarreal
a,
Tomás
González-Lezana
*a and
Paul
Scheier
b
aInstituto de Física Fundamental, IFF-CSIC, Serrano 123, 28006 Madrid, Spain. E-mail: t.gonzalez.lezana@csic.es
bUniversität Innsbruck, Institut für Ionenphysik und Angewandte Physik, Technikerstraße 25, 6020 Innsbruck, Austria
cDipartimento di Chimica, Biologia e Biotecnologie, Universitá di Perugia, 06123 Perugia, Italy
First published on 17th December 2021
Helium clusters around the recently experimentally observed sulphur hexafluoride SF6+ and sulphur pentafluoride SF5+ ions are investigated using a combined experimental and theoretical effort. Mass spectrometry ion yields are obtained and the energetics and structure of the corresponding HeN–SF6+ and HeN–SF5+ clusters are analyzed using path integral molecular dynamics calculations as a function of N, the number of He atoms, employing a new intermolecular potential describing the interaction between the dopant and the surrounding helium. The new force field is optimized on benchmark potential energy ab initio calculations and represented by improved Lennard-Jonnes analytical expressions. This procedure improves the previous potentials employed in similar simulations for neutral SF6 attached to helium nanodroplets. The theoretical analysis explains the characteristic features observed in the experimental ion yields which suggest the existence of stable configurations at specific sizes.
In contrast to all this vast bibliography on sulphur hexafluoride, its ionic counterpart, SF6+, has received considerably less attention. One of the reasons for this is its elusive character which has prevented a reliable experimental detection to date. After a long list of unsuccessful attempts at both observing SF6+15–22 and stabilizing this ionic species in SF6 clusters,23–25 the search for the transient sulphur hexafluoride cation ended with the report of its stabilization using helium nanodroplets (HNDs).26 Albertini et al. demonstrated that sufficiently long-lived SF6+ can be formed by doping charged helium nanodroplets with neutral SF6.26 The ions are identified by means of high-resolution mass spectrometry and collision-induced dissociation following the collision of helium gas with mass-selected HeNSF6+.
Although it seems well accepted that the extra F′ atom of SF6+ (as compared to SF5+) is assumed to occupy an external position forming a weakly bonded complex F′–SF5+,27–30 not too much information regarding the structure of He atoms around the cation can be found in the literature. In their study, Albertini et al. had not paid much attention to the position of the He atoms with respect to SF6+ due to the minor importance of the energetic considerations in elucidating the precise decomposition channels followed in order to produce either the sulphur hexafluoride or pentafluoride cations.26 The authors concede, nevertheless, that the formation of helium snowball cage structures surrounding SF5+ may play a significant role in delaying the recombination of the pair SF5+ + F.26 It is therefore of interest to analyze the pattern followed by He atoms during their solvation of the dopant. In this work, we analyze in detail the clusters formed by He atoms around both cations, SF5+ and SF6+, using path integral molecular dynamics (PIMD) calculations.
For this kind of investigation, the potential describing the interaction between the dopant and the surrounding helium has to be as accurate as possible. However in most of the previous theoretical studies on neutral SF6 attached to HNDs, the intermolecular potential was described using simple empirical two-body (2B) expressions with coefficients and parameters refined by simultaneous fittings of properties such as differential cross-sections, viscosities and virial coefficients.31–34 Some of them were just simple spherical 12-6 or exp-6 potentials with dependence on the experiment considered for the fitting. Typical issues considered in these investigations on the neutral sulphur hexafluoride attached to rare gas droplets are the possible spherical character of octahedral molecules and the manifestation of anisotropic effects.31,35–37 Here, for this work, we have developed a new potential energy surface (PES) to describe the existing molecular interactions. As in previous studies on doped HNDs38–40 benchmark ab initio calculations have been carried out and the resulting potential energy points are optimised and represented via an improved Lennard-Jones (ILJ) expression.
The structure of the paper is as follows. in Section 2 we present the experimental work and in Section 3 we discuss the theoretical part of our study. In particular, Section 3.1 is devoted to the calculation of the intermolecular potential employed in our simulations and in Section 3.3, brief details of the PIMD method are given. Results are shown in Section 4 and conclusions are drawn in Section 5.
The global interaction between the SF5+ ion and an external atom (either F′ or He) can be formulated as a combination of three “effective” components:
Vinter = VvdW + Vind + Velect, | (1) |
The van der Waals VvdW term is expressed as a sum of effective atoms (on SF5+)-external atom pair-wise contributions
(2) |
The formulation adopted for each term Vi term in eqn (2) is of the ILJ type:44
(3) |
(4) |
The key feature of the ILJ functional form is the adoption of additional (variable) n exponential parameters providing more flexibility than the usual Lennard-Jones (12,6) (LJ) ones, thanks to its dependence on Ri as follows:44
n(x) = β + 4.0x2 | (5) |
As for the Vind term of eqn (1), it has been introduced to describe the attractive charge-induced dipole contribution determined by the integer positive charge on SF5+. As suggested by CM545 atomic charge calculations performed at the Hartree–Fock/aug-cc-VTZ46 and B3LYP/cc-VTZ46 levels of theory, here it is assumed that the whole charge is exclusively borne by the S atom and the used expression is
(6) |
In the study reported here, the Velect term of eqn (1) has been introduced only in the case of the SF5+–F′ interaction by retaining the main charge-quadrupole contribution. In particular, the use of the expression
(7) |
In the case of the He–F′ interaction just the first term in eqn (1) is retained since contributions other than the van der Waals one are not necessary and its representation only involves the formula in eqn (3).
All the used parameters are reported in Table 1. Fine tuning has been carried out for ε and Rm by exploiting a comparison with ab initio calculations, performed by using the Molpro2012.1 package,47 of the intermolecular interaction energies. In particular, the optimisation of the force field has been performed by varying the potential parameters within restricted ranges in order to maintain their correct relation with basic properties of involved partners. This guarantees the correctness of the force field represented in the full space of the relative configurations of the interacting partners. For instance, the parameters of the He–F′ pair fall in the right scale of the experimental determination,48 while those of the other pairs scale according to the variation of the electronic polarizability of the interacting partners.49 Results obtained for selected configurations of the interacting partners are shown in Fig. 2–4 for the He–SF5+, F′–SF5+ and He–F′ interactions, respectively. Benchmark theoretical values for the counterpoise corrected interaction energies have been reported and they have been obtained at the CCSD(T) level of theory with two different basis sets (aug-cc-pVQZ and aug-cc-pV5Z46) which allowed the estimation of reliable complete basis set (CBS) extrapolations.50,51 In the calculations involving the SF5+ molecule, the latter has been considered as a rigid body and the used equilibrium geometry is that obtained and reported in the ESI of ref. 26.
F′–SF5+ | ε (meV) | R m (Å) | β | C 4 (meV Å4) | C 3 (meV Å3) |
---|---|---|---|---|---|
F′–S | 4.690 | 3.676 | 8 | −4032 | −2016.2 |
F′–F | 4.540 | 3.272 | 8 |
He–SF5+ | ε | R m | β | C 4 |
---|---|---|---|---|
He–S | 1.914 | 3.556 | 8 | −1440 |
He–F | 2.240 | 3.050 | 8 |
He–F′ | ε | R m | β |
---|---|---|---|
He–F′ | 1.909 | 3.091 | 9 |
Fig. 2 He–SF5+ intermolecular interaction as obtained from ab initio calculations and its analytical representation following eqn (1)–(3). Cuts of the potential (measured in meV) at different values of the angles ϕ and θ as a function of the distance R (in Å) of the He atom with respect to the center of mass of the dopant (see the graphical scheme included in the first panel with the corresponding coordinates). From left to right panels three different positions: (i) hollow 1 with ϕ = 90°, θ = 0°; (ii) hollow 2 with ϕ = 45°, θ = 180° and (iii) top with ϕ = θ = 0°. In the fourth panel, we show the ϕ dependence for R = 3.25 Å (roughly corresponding to the minimum R value for the hollow 1 approach) and, that is, the pathway connecting the hollow 1 and top configurations. |
Fig. 3 Same as Fig. 2 for the F–SF 5+ interaction. |
In general, in Fig. 2–4 a quite good agreement can be observed between the CCSD(T)/CBS ab initio results and the analytical representation of the PESs. More in detail, from Fig. 2 and 3 it can be appreciated that the He-SF5+ and F′-SF5+ interactions show similar features with potential curves providing minima at close intermolecular distances (differences are around 0.1 Å) but with a global interaction being about four times more attractive in the case of the F′ external partner. Moreover, from the last panel of Fig. 2 and 3 it is evident that the most attractive configurations are those with an external atom around the equatorial region (ϕ ∼ 90°) of the SF5+ molecule while locations of the F close to the polar region (ϕ ∼ 0°) provide the least favorable approaches. In fact, the global minimum is found for θ = 0° and ϕ around 65–70 degrees for both He–SF5+ and F′–SF5+ dimers at very close intermolecular distances of around 3.3 Å. Therefore, these results confirm that the F′–SF5+ interaction is of a non-covalent type, even if it is globally quite strong compared to the remaining He–SF5+, He–F′ and He–He contributions.
The present two-body model for the representation of the HeN–SF5+ and HeN–SF6+ interactions can be considered as appropriate, as shown in the ESI† where the negligible role played by three-body effects is analyzed using ab initio calculations. In particular, Fig. S1 of the ESI† reveals a good agreement between the total three-body interaction and that from a pairwise two-body approach as a function of the rotation angle of one He atom with respect to the other one in a He–S–He configuration.
More specifically, initial populations of M individuals or clusters consisting of N He atoms surrounding either the SF5+ or the SF6+ core are generated. Each individual i is characterized by the pair of vectors (i, i) representing the 3N Cartesian coordinates of the atoms and standard deviations for Gaussian mutations, respectively. Initial values of ηi = 1 and random choices for the positions within a specific range (0, Δ), are considered. Thus, each parent creates a single offspring according to
(8) |
(9) |
Pairwise comparisons of the energy of each individual with q random choices as opponents over the union of 2 elements formed of parents (xi, ηi) and offsprings are performed. Individuals with the lowest energy in their competitions with some other opponents are awarded with points, and finally, those individuals out of the union of parents and offsprings with a larger number of winning points are selected as survivors to the next generation, thus becoming new parents. The procedure is repeated until the difference between the potential energies of consecutive generations is lower than the above mentioned tolerance value.
Instead of performing on-the-fly ab initio calculations, a capability included in i-PI, and for saving run time, we used the analytic potential model described above. In this model, the SF5+ core is considered as rigid with the particles arranged at the equilibrium geometry. In order to simulate this behavior in the PIMD runs, the interaction of each pair of particles within the core is described using a very stiff harmonic oscillator with a force constant of 0.01 a.u. This has no consequences for the classical simulation (number of beads M = 1), as the system looks for the minimum of the full PES but needs a separate SF5+ calculation for M > 1 (calibration) to be performed. In fact, for M = 20 (see below), this rigid compound presents a bond energy of 39.03 meV which has to be subtracted from the energies of different HeN-SFn+ complexes.
Based on the white noise Langevin thermostat, we use the global version of the path integral Langevin equation (PILE-G) stochastic thermostatting scheme54 with a unique input parameter τ0, the friction coefficient which determines the strength of the thermostat. For a temperature of 2 K a value of τ0= 1 fs was considered along the simulations. When using a large simulation cubic cell (side = 95 Å) it is not necessary to incorporate barostats as the pressure always remains close to zero. A time interval of Δt = 0.1 fs was chosen to be of the order of 1/500 times, the smallest period in the physical system (∼50 fs, corresponding to the maximum kinetic energy of the SF5+–F′ interaction ∼85 meV), and the quality of the simulation was controlled through the so-called effective energy57 in addition to temperature. The latter oscillates around 2 K within 0.05 K while the former is kept within a variation of ∼0.1%.
We carry out the simulations using optimized minimum energy structures obtained by means of the EA52 as the initial configuration. Details about this method have been given before and, for an example, we invite the interested reader to examine previous applications.53 Thus, starting from those classically estimated minima, we perform the PIMD calculations, first considering a number of beads M in the extended system (ring polymer) M = 1 (classical), and then M = 20. The latter was adopted after using a simple effective atom–atom model for He-SF5+ which leads, by solving the Schrödinger equation, to a binding energy of 14.89 meV (the PIMD value, at M = 20, is 15.19 meV) in such a way that this modest number of beads is able to account for quantum effects, excluding those relative to He–He interactions which would need the use of a huge number of beads (M ≥ 500).
Initial velocities, starting from the initial configurations for the complex produced by EA, were generated from a Maxwell–Boltzmann distribution at a given temperature. All the magnitudes were estimated in the centroid approximation.
Fig. 5 Experimental ion yield for HeNSF5+ for values at gas pressures of 0.16 Pa (blue triangles), 0.18 Pa (red circles) and 0.22 Pa (black squares). |
Minimum energy structures with a He atom with respect to SF5+ according to the present force field global optimization process have been obtained using the EA described in Section 3.2. The optimized configuration for He-SF5+ is shown in panel (a) of the inset of Fig. 4. The classical result reveals the preference of a He atom to occupy a position with values of the spherical angles (see the inset of Fig. 2) of ϕ ∼ 66.5° and θ ∼ 0° and separated at a distance of about R = 3.40 Å from the S atom and about 3.02 Å with respect to the two closest F atoms. There exist, in fact, six symmetrically equivalent minima for either the He atom or the extra F′ atom located at (θ = 0°, ϕ = 66.5°/113.5°), (θ = 120°, ϕ = 66.5°/113.5°), and (θ = 240°, ϕ = 66.5°/113.5°), respectively. We find that this location is also the most favorable site for the extra F′ atom in the sulphur hexafluoride ion.
EA minimization simulations where He atoms are free to move, whereas S and F atoms remain fixed, are carried out to explore the corresponding minimum energy structures for the rest of the HeNSF5+ clusters. Results of the evaporation energies, ΔEN = EN+1 − EN, as a function of the number of He atoms existing in the clusters are shown in Fig. 6. Three noticeable features are seen at specific sizes. In particular, for N = 6, 12 and 24, the ΔEN curve displays a sudden decrease with respect to the almost average value exhibited for immediately smaller clusters. As in previous investigations in doped helium droplets38–40 these effects usually correspond to the filling of a layer or specific caging structures surrounding the dopant. In this case, the analysis of the associated structures for He6SF5+, He12SF5+ and He24SF5+ certainly reveals special arrangements of the He atoms. In particular, for N = 6, all equivalent minimum positions shown in Fig. 4 are occupied. This explains that adding an extra helium atom leads to a decrease of the evaporation energy from ∼20 meV to ∼15 meV as observed in Fig. 6, where an inset of the corresponding He6SF5+ structure is included. Analogously, the addition of six more He atoms yields the construction of an outer cage in which extra helium atoms are located in two triangles with vertices facing the vacant F–Ŝ–F angles. For N = 24, the classical prediction for the minimum energy geometry seems to correspond to a structure formed with 12 He atoms surrounding the above mentioned He12SF5+ as an internal core. The overall appearance could be understood as triangles, both at the top and bottom, and three pairs of He atoms perpendicular to the planar F atoms.
Fig. 6 also confirms that the QM PIMD calculation seems to confirm the presence of the structures predicted by the EA. Despite its more diffuse trend in comparison with the well-defined plateau regions of the classical result, the quantum evaporation energies shown in Fig. 6 show a qualitatively sudden drop for the same sizes.
This apparent success of our present calculations to provide some insights regarding the origin of most of the features observed in the experimental ion yields in Fig. 6, contrasts, nevertheless, with the situation for N = 20, which remains as an intriguing case, since no clear indications of a specific behaviour with respect to consecutive sizes are seen. The analysis of the minimum energy configuration as predicted by EA optimization, shown in Fig. 7, reveals that He atoms keep the minimum configuration observed for N = 12 as an internal core, with the remaining atoms associated with a triangle over the axial straight F–S–F direction, two pairs and an isolated He atom. A perhaps more symmetric closed geometry for such a size of the doped helium cluster is also included in panel (b) of Fig. 7. The N = 12-structure is now surrounded by three pairs of He atoms at the plane formed by the S atom and the three central F atoms and two independent atoms both at the two extremes of the F–S–F axis. Its energy is, however, about 7 meV above the minimum energy geometry shown in panel (a).
Fig. 7 Selected geometries corresponding to He20SF5+. (a) Minimum energy configuration as predicted by EA optimizations and (b) a closed symmetric structure of about 7 meV above the minimum shown in panel (a). As in Fig. 6, He atoms have been artificially bonded as a guide to the eye to distinguish structures around the dopant. |
We have tried to grasp a closer insight by analysing the geometrical structure of He20-SF5+. Fig. 8 shows the energies of this droplet as a function of the PIMD simulation step. We have found that helium atoms form configurations during the simulation which are, in essence, small distortions with respect to precisely the minimum energy configurations obtained with EA optimization. In particular, in Fig. 8, we include the structure obtained with the geometrical centers of the M = 20 beads employed in the calculation for each He atom for that step in which the PIMD energy reaches its minimum. This geometry is basically the same as the optimized classical minimum energy shown above in panel (a) of Fig. 7. This confirms the lack of convergence towards this structure as a possible sufficiently stable configuration. But, aside from this consistency between the quantum mechanical PIMD calculation and the minimum energy geometry found classically with the EA calculation (which, on the other hand, is also seen for the cases of N = 6, 12 and 24) nothing else can be said regarding a special stability for this geometry.
Fig. 9 Same as Fig. 5 for HeNSF6+ with values of gas pressures of 0.11 Pa (green diamonds), 0.16 Pa (blue triangles), 0.18 Pa (red circles) and 0.22 Pa (black squares). |
This suggestion is confirmed when we calculate the evaporation energies by means of the above mentioned EA to search for the location of the He atoms in their minimum energy configurations. In these classical optimizations, the extra F is located in the minimum shown in panel (b) of Fig. 4, whereas the remaining F and S atoms are fixed in their equilibrium locations.
As shown in panel (b) of Fig. 4, the classical optimization applied to one He atom yields the occupancy of the symmetrically equivalent position. Estimates made by Albertini et al. compared different isomers of an atom of He bound to SF6+ (see Fig. 4 from ref. 26). According to the relative energies given in that work, the most stable configuration corresponds to a geometry in which the He atom occupies a different minimum site symmetrically opposed to the extra F atom separated by a central SF5 core (IIa in ref. 26). However, we find that the energy for the He-SF6+ geometry shown in Fig. 4 (IIb in ref. 26) remains about 2.01 meV below that for the isomer IIa. Moreover, the relative energies of the so-called IIc isomer (with the He atom aligned in the axis F–S–F) in that reference would be about 15 meV and 13 meV with respect to the IIb and IIa isomers, respectively (instead of the 3 meV and 11 meV reported in the work by Albertini and coworkers). These results indicate the close proximity of the absolute and relative minima for this system. The slight differences in the actual energy values obtained in both studies are likely to have their origin in the geometry optimization performed in each case. On the one hand the basis set employed in ref. 26 is smaller than the one we use here, and on the other hand, the relaxation of the SF5+ core is not allowed in our approach. Extended calculations by Milan Ončák refining those values published in ref. 26 confirm our present findings.58
The corresponding ΔEN energies as a function of N are shown in Fig. 10. The similarities with the sulphur pentafluoride ion are also manifested in step-like structures suggesting the onset of stable configurations once He atoms fill specific equivalent locations around the dopant. Thus, special features of the evaporation energies are then manifested at N = 5, 11 and 23, that is, exactly at one-He-atom-less sizes as compared with the sulphur pentafluoride ions (see Fig. 6). For those apparent magic numbers, minimum energy geometries as predicted by the EA optimization have been included in the figure. In all these three cases, the extra F atom occupies one of the six equivalent minima found for the He–SF5+ interaction, which, as we mentioned in Section 4.1, was the site reserved for one He atom.
Fig. 10 Same as Fig. 6 for the HeNSF6+ clusters. In this case, the minimum energy structures correspond to N= 5, 11 and 23 droplets. |
The corresponding PIMD simulation yields evaporation energies which, once again, are in a reasonably good qualitative agreement with the classical EA predictions. Although these theoretical results allow therefore the understanding of the presence of the anomalous features at N = 5 and 11 observed in the experimental yields, our calculations do not explain the decrease seen at N = 19 for the different gas pressures.
In an attempt to investigate the location of both the He atoms and the extra F atom in the HeNSF6+ species in more detail, radial probability densities for the different interparticle distances have been obtained. Fig. 11 shows such density functions for the S–F, S–He and He–He distances in the cases of He11–SF6+ and He12–SF6+. The figure includes the comparison between the PIMD results and those of its classical version in which the number of beads M = 1 (see Section 3.3). Radial distributions obtained using Gaussian functions centered at the discrete distances predicted by the EA approach are also compared, thus enabling an overall measurement of the degree of fluctuation of the quantum PIMD result. It is worth noting that in the EA global optimization, the F atom is chosen to remain fixed in the position found for the minimum of the He atoms with respect to the SF5 core. Fig. 11 reveals nevertheless that the quantum distribution for S–F does not differ too much with respect to the classical M = 1 and EA distributions. The same is seen for the S–He distance, with both the PIMD and EA distributions for He12-SF6+ showing that the presence of an extra He atom with respect to N = 11 leads to the onset of a maximum at a slightly larger distance (∼4.4 Å). As expected, more significant discrepancies are observed for the He–He distances, given the weaker interaction between helium atoms as compared to the other components. A similar comparative analysis for some other specific clusters (N = 5–6 and 23–24) is presented in Fig. S2–S4 of the ESI.†
Footnotes |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/d1cp04725f |
‡ Present address: Institute of Food, Nutrition and Health, ETH Zürich, Schmelzbergstrasse 9, 8092 Zürich, Switzerland. |
This journal is © the Owner Societies 2022 |