Very weak interactions: structures, energies and bonding in the tetramers and pentamers of hydrogen sulfide

Abstract

Potential energy surfaces (PESs) for the hydrogen sulfide tetramers and pentamers are shown to be very complex. 11 and 15 different isomers were located on the MP2/6-311++G(d,p) PES of (H₂S)₄ and (H₂S)₅ respectively. CCSD(T) energy calculations on the MP2 optimized geometries suggest that all tetramers are within 2.0 kcal mol⁻¹ of the lowest energy structure, while for the pentamers, all structures are found in a 3.5 kcal mol⁻¹ range. To the best of our knowledge, we report and analyze here for the first time in the scientific literature, a newly found type of very weakly stabilizing intermolecular H₂S⋯SH₂ interaction. In conjunction with traditional H₂S⋯H–S–H hydrogen bonds, these previously unreported H₂S⋯SH₂ intermolecular contacts dictate cluster structures and energies. Our results reveal a very complicated scenario, where a number of different tetramers and pentamers are very close in energy, rendering impossible the unequivocal identification of the global minimum in each case, and as a consequence, suggesting that the properties of these systems would have contributions from many different structures.

---

Introduction

Hydrogen sulfide, H₂S, is a well known poisonous gas of molecular weight 34.08 g mol⁻¹. In nature, hydrogen sulfide is found in a large variety of environments, including volcanic vapors, thermal waters, gas wells, glades, geothermic sources, among others. Commonly, H₂S is produced by the decomposition of organic matter.^1,2 We include next a reduced list of the implications of hydrogen sulfide in a number of important biological, chemical and physical processes:

1. H₂S in endogenously produced in mammals via the transformation of L-cysteine under the catalytic action of cystathionine β-synthase (CBS), cystathionine γ-lyase (CSE),³ and the recently identified 3-mercaptopyruvate sulfurtransferase (3MST).⁴ Deficiencies of the CBS enzyme are known to produce homocystinuria, a pathology associated with mental retardation, sight problems, skeleton abnormalities and episodes of thromboembolism.^5–7 Alongside NO and CO, H₂S is considered to be a physiological gaseous mediator or gaseous transmitter,^8–10 as for example in the vasodilation of the cardiovascular system, facilitating the formation of new blood vessels from pre-existing ones, a process known as angiogenesis.^11–13 H₂S plays a pivotal role in brain functions, acting as neuromodulator and as intracellular messenger.^14–18

2. Hydrogen sulfide is as toxic as cyanide, leading to quick death at concentrations higher than 500 ppm.^19,20 Inhibition of the action of the cytochrome c oxidase enzyme, involved in mitochondrial respiration, has been invoked to explain the lethal action of H₂S leading to cellular death by oxygen deficiency.^20–22 Because of its extreme toxicity, H₂S was used as a chemical weapon during World War I. Very large amounts of H₂S produced because of the metabolism of sulfate-reducing bacteria in ancient oceans are thought to have been responsible for the largest extinction event in the history of this planet. This massive extinction took place during the late Permian period, around 250 million years ago, and it is believed to have killed about 95% of all life on earth, including plants and insects.

3. There is an emerging research line in the modification of drugs aiming at potentiating their pharmacological action. Chattopadhyay and coworkers for example have shown that four H₂S-releasing non steroidal anti inflammatory drugs (H₂S-NSAIDs) inhibit several cancer growth lines in tissues such as pancreas, colon and lung.²³ H₂S-NSAIDs were shown to be 21–3000 times more potent than their NSAIDs counterparts.

4. Hydrogen sulfide has environmental impact. H₂S, being a component of crude oils, is an air pollutant, produced by the burning of petroleum derived fuels;²⁴ the pollution produced by H₂S is a severe hazard, to the extent that desulfuration is a major industrial process in oil refining, in the purification of transportation fuels (gasoline, diesel, aviation fuel), in the preparation of natural gas for domestic use, and so forth.^2,25,26 These processes are important in view of the strict regulations that try to enforce the production of cleaner fuels. SO₂, a sub product of the burning of fuels with high sulfur quantities is one of the main precursors of acid rain.^2,27 In the search for pollutant-free fuels, hydrogen sulfide has been used as a hydrogen source.^1,24,28–31 One of the most used methods for the production of H₂ is the high temperature thermal decomposition of H₂S (ref. 24 and 32) according to the reaction.

Because of a smaller enthalpy of formation (ΔH_f = −4.77 kcal mol⁻¹), H₂S requires 90% less energy than water (ΔH_f = −68.32 kcal mol⁻¹) for the release of H₂.

5. Hydrogen sulfide is a common component of molecular clouds in interstellar space.^33,34 H₂S was the first sulfur compound detected in comets.³⁵ It has also been detected in the volcanic moons of Jupiter, but has not yet been detected in the massive atmosphere of the Jovian planet itself. The search for H₂S in Jupiter is motivated by interest in understanding the composition of its colored clouds. Ammonia polysulfur compounds and other sulfur containing substances are the most promising candidates to explain their puzzling colors.³⁶

Despite the obvious importance in understanding the structures, bonding, and interactions leading to the stabilization of pure H₂S clusters, the scientific literature is completely void of reports concerning clusters beyond (H₂S)₃, with a limited number of either experimental or theoretical works dealing with the dimers and trimers.⁵⁷ This lack of studies of H₂S clusters may be explained by the extreme experimental and theoretical difficulties that arise when treating systems that are held together by very weak interactions. In view of the preceding discussion, in this work, aiming at contributing to the understanding of the structural preferences, energies, and intricate bonding in H₂S clusters, we attempt stochastic explorations and further characterization of the Potential Energy Surfaces of the (H₂S)_n, n = 4, 5 systems. To that end, we locate and characterize local and global minima on both PESs at the MP2/6-311 G++(d,p) level and then use highly correlated CCSD(T)/6-311++G(d,p) calculations to asses energy differences in these highly dispersive systems. In addition, we offer a Bader's Quantum Theory of Atoms in Molecules (QTAIM)³⁷ picture of intermolecular interactions leading to cluster stability.

Computer methods

We used the ASCEC (Annealing Simulado con Energía Cuántica) program^38–40 to produce sets of cluster candidate structures after random explorations of the potential energy surfaces of the (H₂S)_n, n = 4, 5 systems. Details of the algorithm, including descriptions of the production of Markov chains in configuration spaces and the use of a modified Metropolis acceptance test can be found elsewhere.^39,40 The candidate structures were optimized at the MP2/6-311++G(d,p) level without imposition of symmetry constraints; the located stationary points were then characterized as true minima in the corresponding PES by calculating analytical second order derivatives. Single point energy calculations at the CCSD(T)/6-311++G(d,p) level were computed on all MP2 minima. This particular choice of methodology has proven very reliable to treat related hydrogen bonding networks;^39–55 in particular, we expect the CCSD(T) calculations to provide very accurate treatment of the delicate dispersion interactions at play stabilizing the clusters. All optimization, frequency and energy calculations were carried out using the Gaussian09 suite of programs.⁵⁶

Binding energies (BE) were calculated by substracting the energy of a given cluster from n times the energy of the isolated H₂S, in this way, larger positive numbers correspond to larger binding energies. Relative energies (ΔE) were obtained as the difference between the energy of a particular cluster and the energy of the global minimum in the corresponding PES. All MP2, CCSD(T) binding and relative energies reported in this work were corrected for the MP2/6-311++G(d,p) zero-point energies (ZPEs). Isomer populations were estimated via standard Boltzmann distributions. Bonding properties were rationalized using the methods provided by Bader's Quantum Theory of Atoms in Molecules (QTAIM).³⁷ All topological properties of electron densities in this work were calculated using the AIMStudio package.⁵⁸ In this work, we do not correct for basis set superposition error, justification for this practice can be found in the work by David and coworkers⁵³ and in other reports dealing with hydrogen bonded networks.

Results and discussion

ASCEC conditions

In order to ensure thorough explorations of the PESs and that the candidate structures were independent of the initial guess, we ran duplicate ASCEC simulations using the big bang approach, that is, the initial geometry for each ASCEC run was constructed by superimposing all H₂S molecules at the centers of cubic boxes of 7, 8 Å (n = 4) and 8, 9 Å (n = 5) lengths; the systems where then allowed to evolve under the annealing conditions. We used geometrical quenching routes with initial temperatures of 400 K, a constant 5% decrease in temperature, and 100 total temperature steps.

Geometries and structural issues

The stochastic sampling and the proper account for dispersive interactions in this work lead to very complex PESs for the title clusters. At the MP2/6-311++G(d,p) level, we located 11 isomers for the H₂S tetramer and 15 isomers for the pentamer. We use two different notations to describe the clusters: on one hand, a subjective, structure based notation (pyramid, bicycle, etc.) is employed; on the other hand, a simplistic, energy driven notation is also used, in this notation we address the clusters as T_n, P_n with T and P referring to tetramer and pentamer respectively, and the subscript describing the energy ordering of the particular structure; as an illustration of this notation T₁ refers to the most table tetramer, while P₃ refers to the third most stable pentamer. The structures, notations, and relative energies are depicted in Fig. 1 (tetramers) and in Fig. 2 (pentamers). In all cases, intermolecular interactions are always stabilizing, but not strong enough to cause significative changes in the geometries of the constituting monomers. QTAIM shows two types of intermolecular interactions: traditional H₂S⋯H–S–H hydrogen bonds (HBs) and less common H₂S⋯SH₂ long distance contacts, which we call disulfur interactions (DSI).

Fig. 1 Structures, notation and CCSD(T)/6-311++G(d,p)//MP2/6-311++G(d,p) relative energies (kcal mol⁻¹) for the H₂S tetramers. Relative energies are calculated with respect to T₁, the lowest energy structure. QTAIM predicted intermolecular hydrogen and disulfur interactions are represented by dotted lines.

Fig. 2 Structures, notation and CCSD(T)/6-311++G(d,p)//MP2/6-311++G(d,p) relative energies (kcal mol⁻¹) for the H₂S pentamers. Relative energies are calculated with respect to P₁, the lowest energy structure. QTAIM predicted intermolecular hydrogen and disulfur interactions are represented by dotted lines.

Tetramer geometries

Following the suggestions by Acelas and coworkers,⁴¹ the 11 tetramer isomers can be grouped into 4 geometrical motifs defined by the positions of the heavy atoms, regardless of the orientation of hydrogen atoms. The first motif, containing the lowest energy isomer includes the pyramid structures T₁ and T₇, these are the most compact structures; T₁ comprises a HB network in which all molecules act as simultaneous donor/acceptor of hydrogen bonds; T₇, on the other hand exhibits the same geometrical pattern, but contains fewer HBs and three disulfur interactions. The second motif includes T₂, T₅, T₈ and T₉, all bicycle-like structures; among them, T₂, is the second lowest energy structure and T₈ is a quasi heavy atom planar bicycle; the structures in this motif can be clearly separated by their HB/DSI ratios: 4/1, 3/2, 3/2, 1/4 respectively. The third motif contains T₃, T₄ and T₆, all heavy atom planar or quasi planar squares; these are cyclic structures having the same planar network of four S⋯H–S hydrogen bonds, the remaining hydrogen atoms positioning themselves up and down the plane as udud, uudd and uuud respectively; no DSIs are seen in this group. The fourth group contains T₁₀, T₁₁, the highest energy structures, which can be described as heavy atom planar trimers interacting with an additional H₂S molecule; no DSIs are observed in this group either.

Pentamer geometries

The large structural variety in the MP2/6-311++G(d,p) PES for the H₂S pentamers is shown in Fig. 2, were notations and relative energies are also included. As in the case of the tetramers, several distinct geometrical motifs can be realized. The first motif includes P₁, P₂, and P₃, the three lowest energy structures; these correspond to a series of trigonal bipyramids being very close in energy and having the largest numbers of intermolecular contacts (see additional analysis below) with all structures exhibiting at least one DSI in the network of stabilizing interactions. The second motif includes P₄ and P₆, a couple of square based pyramids with 3 and 4 disulfur interactions; interestingly, in P₆, from the perspective of the planar tetramer, the fifth H₂S molecule is attached to the tetramer via four DSIs, all directed towards the same sulfur atom. P₅, P₇ and P₁₀ belong to the next geometrical motif, in which a pyramid of four H₂S molecules interacts with an additional molecule either from the sides (P₅) or from the corners (P₇, P₁₀); P₅ has 3 disulfur interactions while P₇ and P₁₀ have none. The next geometrical motif contains P₈, P₁₁, and P₁₂, a set of quasiplanar tetramers interacting with an extra H₂S molecule; all interactions corresponding to S–H⋯S hydrogen bonds with one single DSI seen in P₁₁. The remaining structures, P₉, P₁₃, P₁₄, P₁₅, characterize each one a geometrical motif on its own, with only 1 DSI in P₉.

Classification of the tetramers and pentamers into geometrical motifs is a subjective exercise which helps rationalizing structural variety among the clusters; it is also somewhat arbitrary in the sense that for larger systems, one structure could belong to more than one motif, and because structures belonging to the same motif are not necessarily closer in energy than structures that are not geometrically related. Our stochastic search of the conformational spaces of (H₂S)_4,5 does not afford linear chains nor ramifications with more than one extra molecule on the sides; this is a testament of the weakness of the interactions under study, which lead to a preference of compact structures over open structures with fewer intermolecular contacts. It is also important to notice that the present study should not be considered as definitive for complete characterization for both PESs, this is not an attainable goal with present day experimental nor theoretical methods, due to the overwhelming number of structural possibilities that increase exponentially with the size of the systems; many of the structures reported here for example would have large numbers of virtually isoenergetic isomers because of the relative positions of hydrogen atoms not taking part in the stabilizing HB networks; in addition, the structures reported here are local minima in the corresponding PES within the harmonic approximation for vibrational frequencies, it is to be expected that for very weak interactions, anharmonicities may play an important role in the precise characterization of the intermolecular interactions, however, those calculations are beyond the scope of this work.

Radial distributions

The weak interactions in the H₂S tetramers and pentamers lead to very long hydrogen bonds; we emphasize that the intermolecular contacts reported here are not selected upon visual inspection of the structures in Fig. 1 and in Fig. 2, instead, they correspond to those bonding interactions derived from topological analysis of the electron densities under QTAIM. Fig. 3 shows a distribution of H₂S⋯H–S–H and H₂S⋯SH₂ distances. Hydrogen bonding in the H₂S tetramers and pentamers covers the 2.6–3.6 Å range with a strong peak centered around ≈2.7 Å; the 1.0 Å range for hydrogen bonding leads to a rich structural variety within the same type of interactions, in contrast, the stronger interactions responsible for hydrogen bonding in typical water clusters produce a peak centered around ≈1.9 Å.⁵⁴ Noticeably, the distribution plots for tetramers and pentamers have the same shapes, differing only in the number of interactions. Disulfur interactions cover the 3.6–4.2 Å interval, with a comparatively reduced, but not negligible number of occurrences.

Fig. 3 Radial distributions (Å) for H₂S⋯HSH and H₂S⋯SH₂ interactions in the MP2/6-311++G(d,p) potential energy surfaces of the H₂S tetramers and pentamers. Pointed and dashed lines correspond to hydrogen bond distances in the tetramers and pentamers respectively, the solid line reflects the total number of hydrogen bonds occurring at a given distance regardless of the size of the cluster.

Energies, cluster stabilization and bonding

Table 1 summarizes the most relevant energy results calculated for the tetramers and pentamers of hydrogen sulfide: binding and relative energies at the MP2/6-311++G(d,p) and CCSD(T)/6-311++G(d,p)//MP2/6-311++G(d,p) levels are listed in conjunction with estimated Boltzmann populations for all structures reported here. Two very important observations are realized upon inspection of Table 1: first, the interactions holding the clusters as discrete units are very weak, producing binding energies no larger than 6.88 and 9.80 kcal mol⁻¹ for the tetramers and pentamers respectively, in contrast, replacing the sulfur atom by an oxygen atom results in binding energies for the water tetramers and pentamers calculated at the same levels of theory that exceed 30 kcal mol⁻¹;^39,43 second, structural isomers are separated by very small energy gaps, leading to a complex situation where all structures have significant populations. Since as discussed above, each structure allows a large number of almost isoenergetic isomers because of the relative positions of hydrogen atoms, it follows that properties of H₂S tetramers and pentamers would have contributions from multitudes of isomers; under these conditions, unequivocal identification of the global minima is not possible.

Table 1 Energetic analysis (kcal mol⁻¹) for the (H₂S)_4,5 clusters. ΔE: relative energy with respect to the lowest energy structure. BE : ZPE corrected binding energy^a

Structure	ΔE^b	ΔE^c	BE^b	BE^c	%xi^b	%xi^c	%xi^d
a %xi: estimated Boltzmann populations.b MP2/6-311++G(d,p).c CCSD(T)/6-311++G(d,p)//MP2/6-311++G(d,p).d MP2/6-311++G(d,p) Gibbs free energies (298 K, 1 atm).
T1, pyramid a	0.00	0.00	6.88	6.35	42.0	38.8	13.4
T2, bicycle a	0.83	0.77	6.05	5.59	10.2	10.7	3.3
T3, planar square a	0.79	0.88	6.09	5.47	10.2	9.4	3.3
T4, planar square b	0.94	1.00	5.94	5.35	9.0	7.3	6.2
T5, bicycle b	1.17	1.00	5.71	5.35	6.1	7.3	5.4
T6, planar square c	1.00	1.06	5.88	5.29	7.9	7.3	5.4
T7, pyramid b	1.30	1.09	5.58	5.26	4.7	6.4	0.6
T8, bicycle c	1.41	1.25	5.47	5.10	3.7	5.0	28.9
T9, bicycle d	1.73	1.46	5.15	4.89	2.2	3.4	1.0
T10, trimer + 1a	1.76	1.69	5.11	4.66	2.2	2.3	17.3
T11, trimer + 1b	1.88	1.81	5.06	4.54	1.7	2.0	15.2
P1, trigonal bipyramid a	0.00	0.00	9.80	9.24	21.7	23.4	10.9
P2, trigonal bipyramid b	0.07	0.10	9.73	9.14	19.1	20.6	3.0
P3, trigonal bipyramid c	0.26	0.31	9.54	8.93	14.7	14.0	0.8
P4, square pyramid a	0.25	0.35	9.55	8.89	14.7	14.0	3.9
P5, pyramid + 1a	0.51	0.60	9.30	8.64	8.8	8.4	1.4
P6, square pyramid b	0.53	0.64	9.28	8.60	8.8	8.4	0.4
P7, pyramid + 1b	1.41	1.44	8.40	7.80	2.4	2.0	5.0
P8, tetramer/trimer a	1.29	1.48	8.52	7.76	2.1	2.0	10.9
P9, tricycle	1.54	1.53	8.26	7.71	1.7	1.8	8.4
P10, pyramid + 1c	1.48	1.55	8.32	7.69	1.9	1.8	18.2
P11, tetramer/trimer b	1.70	1.85	8.10	7.39	1.3	1.1	1.2
P12, tetramer/trimer c	1.86	1.92	7.94	7.32	1.0	0.9	0.8
P13, double trimer	1.92	2.03	7.88	7.21	0.9	0.8	3.0
P14, puckered ring	1.99	2.20	7.81	7.04	0.8	0.6	5.0
P15, trimer + 1 + 1 trans	3.45	3.43	6.35	5.81	0.1	0.1	25.8

---

It has been reported that MP2 overestimates the strength of interactions in dispersive complexes.^59,60 This may have to do with the limitations imposed to truncated perturbative expansions when the perturbation, as in the case of a considerably large number of interacting electrons, is not too small. In our case, Table 1 shows an excellent agreement between MP2 and CCSD(T) energies: MP2 binding energies for example are not overestimated by more than 0.77 kcal mol⁻¹ (roughly ≈10% of the total BE).

On the grounds of purely ZPE corrected electronic energies, Table 1, Fig. 1, and 2 suggest a preference for compact, three dimensional structures for both H₂S tetramers and pentamers. This is in clear contrast with structural preferences for water clusters of similar molecularities, where three dimensional structures appear to be favored only for n ≥ 6;^39,41,43,54 such structural preferences are related to the presence of cage critical points (CCPs) in the electron densities only for the hexamer,⁵⁵ and larger water clusters; in our case, T₁ and P₁, the lowest energy tetramer and pentamer exhibit at least one CCP. On the other hand, Gibbs energies calculated at 298 K, 1 atm, (Table 1) reveal that once temperature, entropy, and internal degrees of freedom are considered, more open structures are preferred; this situation has already been pointed out for other types of small clusters^41,54 and suggest that at high enough temperatures, the systems approach the ideal gas, where intermolecular interactions are not favored.

Topological analysis of the electron densities

In order to study the nature of chemical bonding in the weakly bound tetramers and pentamers of H₂S, we used the Quantum Theory of Atoms In Molecules (QTAIM) developed by Bader;³⁷ in particular, we gain insight from the properties of critical points in the electron densities. Critical points are characterized by vanishing gradients →∇ρ(r_c) = 0; besides nuclear and non-nuclear attractors, such points are commonly classified in three categories: bond (BCPs), ring (RCPs) and cage (CCPs) critical points.

QTAIM suggests that H₂S⋯H–S–H hydrogen bonds (HBs) and H₂S⋯SH₂ disulfur interactions (DSIs) are the only interactions responsible for cluster stabilization. Therefore, within the QTAIM framework, in this paper, we describe bonding in the (H₂S)_4,5 clusters from the properties obtained at hydrogen and disulfur bond critical points; such properties include the electron density ρ(r_c) and its Laplacian ∇²ρ(r_c), as well as the potential [scr V, script letter V](r_c), kinetic [capital G, script](r_c), and total [script letter H](r_c) energy densities. QTAIM describes bonding as the result of contributions from closed shell, intermediate and shared interactions; in this work, we quantify the relative strength of HBs and DSIs using the criteria proposed by Espinosa and coworkers⁶¹ summarized in the following scheme:

Table 2 shows that larger numbers of HBs or larger numbers of DSIs do not necessarily imply larger interaction energies: it is seen that cluster stabilization is qualitatively related, but is not directly correlated with the number of intermolecular interactions nor with the topological complexity of the electron densities; take for example the two pyramid structures, T₁ and T₇ (Fig. 1 and Table 2), they have the same topological complexity and the same number of stabilizing interactions, however, they are separated by more than 1 kcal mol⁻¹ and by five other structures with different geometries and topologies in between, with minute energy differences. Similar lack of correlations between the number of stabilizing interactions and stabilization energies have been reported for the related water tetramers³⁹ and pentamers,⁴³ however, in the present case, the situation is exacerbated because of the very small energy differences among isomers.

Table 2 Topology of the electron densities for the hydrogen sulfide tetramers and pentamers. HBCP, RCP, CCP and DSICP correspond to the number of hydrogen bond, ring, cage, and disulfur interaction critical points respectively, ΣCP is the sum of all of them. ΔE are the CCSD(T)//MP2 relative energies taken from Table 1

Structure	ΔE	HBCP	RCP	CCP	DSICP	ΣCP
T1, pyramid a	0.00	6	4	1	0	11
T2, bicycle a	0.77	4	2	0	1	7
T3, square planar a	0.88	4	1	0	0	5
T4, square planar b	1.00	4	1	0	0	5
T5, bicycle b	1.00	3	2	0	2	7
T6, square planar c	1.06	4	1	0	0	5
T7, pyramid b	1.09	2	4	1	4	11
T8, bicycle c	1.25	3	2	0	2	7
T9, bicycle d	1.46	1	2	0	4	7
T10, trimer + 1 + 1a	1.69	4	1	0	0	5
T11, trimer + 1 + 1b	1.81	4	1	0	0	5
P1, trigonal bipyramid a	0.00	7	7	2	2	18
P2, trigonal bipyramid b	0.10	7	7	2	2	18
P3, trigonal bipyramid c	0.31	8	7	2	1	18
P4, square pyramid a	0.35	5	4	0	3	12
P5, pyramid + 1a	0.60	5	5	1	3	14
P6, square pyramid b	0.64	4	4	0	4	12
P7, pyramid + 1b	1.44	7	4	1	0	12
P8, tetramer/trimer a	1.48	6	2	0	0	8
P9, tricycle	1.53	6	3	0	1	10
P10, pyramid + 1c	1.55	6	4	1	1	12
P11, tetramer/trimer b	1.85	5	2	0	1	8
P12, tetramer/trimer c	1.92	6	2	0	0	8
P13, double trimer	2.03	6	2	0	0	8
P14, puckered ring	2.20	5	1	0	0	6
P15, trimer + 1 + 1 trans	3.43	5	1	0	0	6

---

In the study of several S–S formal bonds, Knop and coworkers⁶² found a correlation between the electron densities at the bond critical points with S–S bond lengths, their findings led them to state “This relationship, the analogue of which has been demonstrated for Al–F, O–O, and Be–Cl bonds and which is expected to be of general validity, makes possible estimates of r_e, from ρ(r_c) and, conversely, estimates of ρ(r_c), bond order, and related properties from r_e”. As noted by the authors, similar relationships were previously shown to exist for other types of bonds, thus, they postulated that [r_e;ρ(r_c)] relationships (equilibrium bond length, electron density at the bond critical point) exist for any class of binary bonds, regardless of the origin and magnitude of the individual contributions from particular molecular orbitals to the electron densities at BCPs. A later study considering many types of bonds⁶³ supports this hypothesis, reporting logarithmic [r_e;ρ(r_c)] correlations with larger slopes for hydrogen bonds than for covalent bonds. Results obtained studying a large set of systems exhibiting N–H⋯N hydrogen bonds^64,65 also agree with these findings: additional analysis of the N–H⋯N interactions uncovered among many other pairs, strong [r_e;ρ(r_c)] and [r_e;∇²ρ(r_c)] correlations, and suggested that covalent and hydrogen bonds are better described by separate fittings. In this work, to the best of our knowledge, we present for the first time (Fig. 4, top right panel) evidence that intermolecular long distance S⋯S interactions abide to the same type of logarithmic [r_e;ρ(r_c)] relationship, these decreasing logarithmic tendencies accurately describe the asymptotic behavior for large atom separations that eventually leads to no interactions at sufficiently large distances (earlier estimates of linear relationships for other systems would eventually lead to negative densities at BCPs for sufficiently long interaction distances). Fig. 4 clearly shows that HBs and DSIs have the same type of [r_e;ρ(r_c)] correlation but can be distinctively characterized as two physically different interactions: fitting the data to a straight line for HBCPs afforded ln ρ(r_c) = −1.7050r + 0.0828 (R² = 0.98) while for DSICPs the obtained adjusted line is ln ρ(r_c) = −1.3063r + 0.0817 (R² = 0.90). Bond critical points with larger electron densities imply that the involved atoms share electrons to a larger extent in the region surrounding the critical points, leading to larger covalent character for the interactions; accordingly, the top right panel of Fig. 4 shows that electron densities are larger for hydrogen bonding than for the S⋯S interactions, indicating the highly closed shell nature of the intermolecular long distance disulfur contacts.

Fig. 4 Plots of selected quantities calculated at the critical points of the electron densities for the H₂S tetramers and pentamers. (Top left) Typical structure showing that only those QTAIM determined intermolecular interactions are taken into account: solid arrows point to hydrogen bonds, open arrows point to disulfur interactions. (Top right) Logarithmic relationship between bond lengths and electron densities at BCPs. (Bottom left) Quantification of interaction strengths using Espinosa's criteria (eqn (2)). (Bottom right) Potential, kinetic and total energy densities as a function of bond lengths. Data taken from the MP2/6-311++G(d,p) PESs.

Espinosa and coworkers⁶¹ proposed a method for quantification of the covalency character of chemical interactions by local application of the virial theorem to BCPs (eqn (2)). Their model divides the allowed values for |[scr V, script letter V]|/[capital G, script], the ratio of potential to kinetic energy densities evaluated at bond critical points, into three distinct regions: for values in the [0, 1] interval, interactions are to be characterized as closed shell (very little electron sharing), values larger than two correspond to increasingly covalent interactions, while the [1, 2] interval comprises interactions of intermediate character. The bottom left panel of Fig. 4 plots the |[scr V, script letter V]|/[capital G, script] ratio against the bond parameter for the H₂S tetramers and pentamers. Based on this plot, we can state a quantified order of covalency (or of closed shell nature if preferred) consistent with the observations drawn from the top right panel of the same figure pointed out above for all interactions resulting in intermolecular BCPs:HBs are more covalent than DSIs. This quantification clearly places all interactions below the covalency threshold, with S⋯S interactions having a stronger closed shell nature.

The bottom right panel in Fig. 4 plots the potential, kinetic and total energy densities at the BCPs for all types of interactions as a function of the separation between atoms. It is clearly seen that kinetic energy is always a repulsive term, while energetic stabilization to the critical points is provided by the potential energy. A very interesting observation is that total energy densities at all HBCPs and DSICPs are positive, this has been reported in other hydrogen bonded clusters^41,43,46,47 but is not observed in other systems stabilized by closed shell interactions, specifically in small gold clusters.⁶⁶

The H₂S⋯SH₂ interactions

Distances for H₂S⋯SH₂ contacts are quite large, they can be found in the 3.6–4.2 Å interval (Fig. 3). Similarly, S⋯H distances in hydrogen bonds in H₂S tetramers and pentamers are uncommonly large (2.6–3.6 Å Fig. 3), thus, the title clusters exhibit very small interaction energies. To the best of our knowledge, H₂S⋯SH₂ intermolecular interactions are reported for the first time in this work. Our results suggest that disulfur contacts comprise a novel type of very weak stabilizing interactions, we support this claim based in the following evidence:

1. HB/DSI, the ratio of hydrogen bonds to disulfur interactions is smaller than 1 in some clusters (1/4 in T₉ and 2/4 in T₇, Table 2). If DSIs were not stabilizing interactions, the collective destabilization of 4 DSIs would work against the 1 and 2 weakly stabilizing hydrogen bonds so that T₉ and in T₇ would not be local minima in the PES.

2. In P₆ (Fig. 2), the molecule located above the planar tetramer is attached to all four molecules in the plane by disulfur interactions, without a single hydrogen bond; clearly, the only reason why this particular structure is stable is because disulfur interactions are stabilizing.

3. If DSIs were destabilizing, they would be present in larger numbers in the least stable structures, however, the opposite is seen in Table 2: DSIs are encountered among the most stable structures and are not present at all in the least stable tetramers and pentamers.

4. Cluster stability seems to be related to the total number of intermolecular contacts, regardless of their nature.

5. HBs and DSIs have the same type of [r_e;ρ(r_c)] correlation (top right panel in Fig. 4), thus can be distinctively characterized as two physically different stabilizing interactions.

6. DSIs are characterized as closed shell weakly stabilizing interactions (bottom left panel in Fig. 4).

Summary, conclusions and perspectives

Stochastic explorations of the potential energy surfaces for the hydrogen sulfide tetramers and pentamers afforded 11 and 15 isomers respectively at the MP2/6-311++G(d,p) level. Highly correlated CCSD(T) energy calculations on the MP2 optimized structures reveal very small binding energies, not exceeding 7 kcal mol⁻¹ in the tetramers and 10 kcal mol⁻¹ in the pentamers; within these energy ranges, isomers are very close in energy, separated by no more than 2.0 and 3.5 kcal mol⁻¹ respectively, making it impossible to unambiguously identificate the global minima on each PES, this situation leads to contributions from several structures to the properties of H₂S tetramers and pentamers. On the basis of purely ZPE corrected electronic energies, three dimensional, compact structures are energetically preferred, however, consideration of temperature, entropy and internal degrees of freedom when calculating Gibbs free energies at 298 K, 1 atm, suggest a preference for more open structures; this result is consistent with the ideal gas situation, where intermolecular interactions are disfavored at sufficiently high temperatures. QTAIM predicts two types of physically distinct interactions: traditional H₂S⋯H–S–H hydrogen bonds and, to the best of our knowledge, a new type of stabilizing intermolecular H₂S⋯SH₂ contacts which we address as disulfur interactions. Hydrogen bonding in the title systems covers the 2.6–3.6 Å range for S⋯H distances, which is considerably larger than for other hydrogen bonded systems and constitute another aspect of the weakness of the particular interactions treated in this study. Disulfur interactions, on the other hand, spread along the 3.6–4.2 Å. Various criteria derived from topological analyses of the electron densities under the QTAIM framework place all intermolecular interactions in the closed shell category, suggesting that disulfur interactions have a larger degree of “closed shellness” than H₂S⋯H–S–H hydrogen bonds. To the best of our knowledge, we present here for the first time, evidence that suggests that intermolecular long distance S⋯S bonds show logarithmic relationships between the electron density at critical points with the distance separating the atoms.

The findings in this work should not be taken as definitive for complete characterizations of the potential energy surfaces in the hydrogen sulfide tetramers and pentamers, this goal is not realizable with present days experimental nor theoretical techniques due to the very large number of isomers to be found because of the relative positions of hydrogen atoms not involved in the network of stabilizing interactions, nonetheless, we provide a comprehensive spectrum of possibilities and are confident that any missing isomers will fall into the geometrical motifs described within the text.

Acknowledgements

This work was partially supported by Universidad de Antioquia, via “Estrategia de Sostenibilidad 2013–2014”. Very useful comments to the manuscript by Prof. Frank Weinhold, University of Wisconsin, are acknowledged.

References

1. S. Ouali, S. Chader, M. Belhamel and M. Benziada, Int. J. Hydrogen Energy, 2011, 36(6), 4103–4109.
2. C. Sung, S. Al Hashimi, A. McCormick, M. Tsapatsis and M. Cococcioni, J. Phys. Chem. C, 2012, 116(5), 3561–3575.
3. M. Stipanuk and P. Beck, Biochem. J., 1982, 206(2), 267–277.
4. N. Shibuya, M. Tanaka, M. Yoshida, Y. Ogasawara, T. Togawa, K. Ishii and H. Kimura, Antioxid. Redox Signaling, 2009, 11(4), 703–714.
5. B. Fowler, J. Kraus, S. Packman and L. E. Rosenberg, J. Clin. Invest., 1978, 61(3), 645.
6. S. H. Mudd, F. Skovby, H. L. Levy, K. D. Pettigrew, B. Wilcken, R. E. Pyeritz, G. Andria, G. H. Boers, I. L. Bromberg and R. Cerone, Am. J. Hum. Genet., 1985, 37(1), 1.
7. K. Qu, S. Lee, J. Bian, C. Low and P. Wong, Neurochem. Int., 2008, 52(1), 155–165.
8. R. Palmer, D. S. Ashton and S. Moncada, Nature, 1988, 333(6174), 664–666.
9. T. Dawson, V. Dawson and S. H. Snyder, Ann. Neurol., 1992, 32(3), 297–311.
10. L. Li, A. Hsu and P. K. Moore, Pharmacol. Ther., 2009, 123(3), 386–400.
11. H. Kimura, Mol. Neurobiol., 2002, 26(1), 13–19.
12. J. Koenitzer, T. Isbell, H. Patel, G. Benavides, D. Dickinson, R. Patel, V. Darley-Usmar, J. Lancaster, J. Doeller and D. Kraus, Am. J. Physiol.: Heart Circ. Physiol., 2007, 292(4), 1953–1960.
13. L. Rochette and C. Vergely, Ann. Cardiol. Angeiol., 2008, 57(3), 136–138.
14. P. Moore, M. Bhatia and S. Moochhala, Trends Pharmacol. Sci., 2003, 24(12), 609–611.
15. Y. Kimura, Y. Mikami, K. Osumi, M. Tsugane, J. Oka and H. Kimura, FASEB J., 2013, 27(6), 2451–2457.
16. H. Kimura, Exp. Physiol., 2011, 96(9), 833–838.
17. B. H. Tan, P. T. Wong and J. S. Bian, Neurochem. Int., 2010, 56(1), 3–10.
18. E. Lowicka and J. Beltowski, Pharmacol. Rep., 2006, 59(1), 4–24.
19. D. Truong, M. Eghbal, W. Hindmarsh, S. Roth and P. J. O'Brien, Drug Metab. Rev., 2006, 38(4), 733–744.
20. R. Reiffenstein, W. Hulbert and S. Roth, Annu. Rev. Pharmacol. Toxicol., 1992, 32, 109–134.
21. A. Khan, M. Schuler, M. Prior, S. Yong, R. Coppock, L. Florence and L. Lillie, Toxicol. Appl. Pharmacol., 1990, 103(3), 482–490.
22. C. Cooper and G. Brown, J. Bioenerg. Biomembr., 2008, 40(5), 533–539.
23. M. Chattopadhyay, R. Kodela, N. Nath, Y. Dastagirzada, C. Velazquez-Martinez, D. Boring and K. Kashfi, Biochem. Pharmacol., 2012, 83(6), 715–722.
24. S. Baykara, E. Figen, A. Kale and T. Nejat Veziroglu, Int. J. Hydrogen Energy, 2007, 32(9), 1246–1250.
25. R. Yang, A. Hernández-Maldonado and F. Yang, Science, 2003, 301(5629), 79–81.
26. A. Hernández-Maldonado and R. Yang, Catal. Rev.: Sci. Eng., 2004, 46(2), 111–150.
27. H. Song, M. Park, S. Kwon, K. Yi, E. Croiset, Z. Chen and S. Nam, Appl. Surf. Sci., 2013, 276, 646–652.
28. E. Fletcher, J. Noring and J. Murray, Int. J. Hydrogen Energy, 1984, 9(7), 587–593.
29. T. Kappauf and E. Fletcher, Energy, 1989, 14(8), 443–449.
30. C. Linkous, T. Mingo and N. Muradov, Int. J. Hydrogen Energy, 1994, 19(3), 203–208.
31. J. Zaman and A. Chakma, Int. J. Hydrogen Energy, 1995, 20(1), 21–28.
32. Y. Lu and L. Schaefer, J. Power Sources, 2004, 135(1), 184–191.
33. N. Ukita and M. Morris, Astron. Astrophys., 1983, 121, 15–18.
34. P. Thaddeus, M. Kutner, A. Penzias, R. Wilson and K. Jefferts, APJ, 1972, 176, L73.
35. J. Crovisier, D. Despois, D. Bockelée-Morvan, P. Colom and G. Paubert, Icarus, 1991, 93(2), 246–258.
36. J. Joiner, P. G. Steffes and K. S. Noll, IEEE Trans. Microwave Theory Tech., 1992, 40(6), 1101–1109.
37. R. Bader, Atoms in Molecules. A Quantum Theory, Oxford University Press, NY, 1990.
38. J. Pérez and A. Restrepo, ASCEC V–02: Annealing Simulado con Energía Cuántica. Property, development and implementation: Grupo de Química–Física Teórica, Instituto de Química, Universidad de Antioquia, Medellín, Colombia, 2008.
39. J. Pérez, C. Hadad and A. Restrepo, Int. J. Quantum Chem., 2008, 108(10), 1653–1659.
40. J. Pérez, E. Flórez, C. Hadad, P. Fuentealba and A. Restrepo, J. Phys. Chem. A, 2008, 112(25), 5749–5755.
41. N. Acelas, G. Hincapié, D. Guerra and A. Restrepo, J. Chem. Phys., 2013, 139, 044310.
42. J. Murillo, J. David and A. Restrepo, Phys. Chem. Chem. Phys., 2010, 12(36), 10963–10970.
43. F. Ramírez, C. Hadad, D. Guerra, J. David and A. Restrepo, Chem. Phys. Lett., 2011, 507, 229–233.
44. C. Ibargüen, M. Manrique-Moreno, D. Guerra, J. David and A. Restrepo, Phys. Chem. Chem. Phys., 2013, 15, 3203–3211.
45. S. Gómez, D. Guerra, J. David and A. Restrepo, J. Mol. Model., 2013, 19, 2173–2181.
46. D. Guerra, J. David and A. Restrepo, Phys. Chem. Chem. Phys., 2012, 14, 14892.
47. D. Guerra, J. David and A. J. Restrepo, J. Comput. Methods Sci. Eng., 2014, 14, 93.
48. J. Gonzalez, E. Florez, J. Romero, A. Reyes and A. Restrepo, J. Mol. Model., 2013, 19, 1763.
49. S. Mejía, J. Espinal, A. Restrepo and F. Mondragón, J. Phys. Chem. A, 2007, 111, 8250.
50. A. Zapata-Escobar, M. Manrique-Moreno, D. Guerra, C. Hadad and A. Restrepo, J. Chem. Phys., 2014, 140, 184312.
51. J. Romero, A. Reyes, J. David and A. Restrepo, Phys. Chem. Chem. Phys., 2011, 13, 15264.
52. C. Hadad, A. Restrepo, S. Jenkins, F. Ramirez and J. David, Theor. Chem. Acc., 2013, 132, 1376.
53. J. David, D. Guerra and A. Restrepo, J. Phys. Chem. A, 2009, 113, 10167.
54. G. Hincapié, N. Acelas, M. Castaño, J. David and A. Restrepo, J. Phys. Chem. A, 2010, 114, 7809.
55. S. Jenkins, A. Restrepo, J. David, D. Yin and S. Kirk, Phys. Chem. Chem. Phys., 2011, 13, 11644.
56. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. Montgomery Jr, J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, Ö. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski, and D. J. Fox, Gaussian 09, Revision D.01, Gaussian, Inc., Wallingford CT, 2009.
57. H. Tsujii, K. Takizawa and S. Koda, Chem. Phys., 2002, 285, 319.
58. T. Keith, AIMALL version 10.09.12, 2010, http://aim.tkgristmill.com/.
59. S. Tsuzuki and T. Uchimaru, Curr. Org. Chem., 2006, 107, 745–762.
60. T. Tauer, M. Derrick and D. Sherrill, J. Phys. Chem. A, 2005, 109, 191–196.
61. E. Espinosa, I. Alkorta, J. Elguero and E. Mollins, J. Chem. Phys., 2002, 117, 5529.
62. O. Knop, R. Boyd and S. Chos, J. Am. Chem. Soc., 1988, 110, 7299–7301.
63. I. Alkorta, I. Rozas and J. Elguero, Struct. Chem., 1998, 9, 243.
64. O. Knop, K. Rankin and R. Boyd, J. Phys. Chem. A, 2001, 105, 6552.
65. O. Knop, K. Rankin and R. Boyd, J. Phys. Chem. A, 2003, 107, 272.
66. J. David, D. Guerra and A. Restrepo, Chem. Phys. Lett., 2012, 539, 64.

---

Footnote

† Electronic supplementary information (ESI) available: Cartesian coordinates for all tetramers and pentamers reported in this work. See DOI: 10.1039/c4ra09430a

---