Alba Vargas-Caamala,
Sudip Panb,
Filiberto Ortiz-Chic,
Jose Luis Cabellosa,
Roberto A. Botod,
Julia Contreras-Garciad,
Albeiro Restrepo*e,
Pratim K. Chattaraj*b and
Gabriel Merino*a
aDepartamento de Física Aplicada, Centro de Investigación y de Estudios Avanzados, Unidad Mérida. km 6 Antigua Carretera a Progreso. Apdo. Postal 73, Cordemex, 97310, Mérida, Yuc., Mexico. E-mail: gmerino@mda.cinvestav.mx
bDepartment of Chemistry and Center for Theoretical Studies, Indian Institute of Technology, Kharagpur, 721302, India. E-mail: pkc@chem.iitkgp.ernet.in
cInstituto Tecnólogico Superior de Calkiní, Av. Ah-Canul s/n, Carr. Fed. Calkiní-Campeche, CP 24900, Calkiní, Campeche, Mexico
dLaboratoire de Chimie Théorique and CNRS, UMR 7616, 75252 Paris Cedex 05, Paris, France
eInstituto de Química, Universidad de Antioquia UdeA, Calle 70 No. 52-21, Medellín, Colombia. E-mail: albeiro.restrepo@udea.edu.co
First published on 11th November 2015
An exhaustive exploration of the potential energy surfaces of ferrocene, ruthenocene and osmocene dimers has been performed. Our computations involving dispersion show that only four different isomers are present in each metallocene dimer. The collective action of small interaction energies of dispersive nature leads to a dissociation energy of 7.5 kcal mol−1 for the ferrocene dimer. Dispersion has strong effects on the geometrical parameters, reducing the M⋯M distances by almost 1 Å. Our results also reveal that inclusion of entropic factors modifies the relative stability of the complexes. The nature of bonding is examined using the energy decomposition analysis and the non-covalent interaction index. Both analyses indicate that dispersion is the major contributing factor in stabilizing a metallocene dimer.
Recently, Bogdanović and Novaković6 evaluated the frequency of occurrence of the ferrocene dimer in crystals reported in the Cambridge Structural Databank.7 They found that 46.8% of ferrocene derivative crystals contain a dimer where both units are in a parallel orientation, with one of the ferrocene units shifted along the z-axis by half of the Cp–Fe bond length (see Fig. 1).
Fig. 1 Building block identified by Bogdanović and Novaković in ferrocene derivative crystals.6 |
Keeping in mind the study of Bogdanović and Novaković, a detailed analysis of the bonding and nature of the ferrocene dimer will be worth pursuing in order to understand the pattern of its supramolecular arrangements.8 In the present work, we exhaustively explored the potential energy surfaces (PESs) of the ferrocene dimers. In addition, we extended this search to the PESs of ruthenocene and osmocene dimers. Our Kohn–Sham density functional theory (KS-DFT) computations involving dispersion show that only four different isomers are present in each metallocene dimer in the gas phase. Interestingly, the collective action of small interaction energies of dispersive nature leads to a higher dissociation energy for the ferrocene dimer (7.5 kcal mol−1) at low temperatures than that computed for a typical hydrogen bond in the water dimer (4.8 kcal mol−1). Our calculations also reveal that inclusion of entropic factors modifies the relative stability of the complexes. The nature of bonding is examined using the energy decomposition analysis (EDA)9 and the non-covalent interaction (NCI) index.10 Particularly, EDA shows a delicate balance between dispersion and electrostatic contributions.
Our computational procedure employs a modified-kick heuristic algorithm implemented in Bilatu to systematically explore the potential energy surfaces (PESs) of molecular clusters.12 The readers interested in details about this metaheuristic are referred to ref. 13. Final equilibrium geometries are reported at the PBE14/def2-TZVP15 and PBE-D2/def2-TZVP levels. To describe the scalar relativistic effects, an effective core potential was used for Ru and Os, describing the behavior of 28 and 60 core electrons, respectively. The latter approach includes the D2 version of Grimme’s dispersion corrections.16 The energy differences discussed here include the harmonic zero point energy (ZPE) correction. In order to evaluate the thermal effects, we use the procedure described by Irikura17 as is implemented in thermo code,18 where the standard molar entropy and enthalpy changes are computed from the molecular partition function. All the quantities needed are taken from the harmonic vibrational frequency computations. All the computations are done using the Gaussian 09 program.19
The nature of the interactions is analyzed by EDA20 at the revPBE-D321/TZ2P//PBE-D2/def2-TZVP level using the ADF (2013.01) package.22 We do not use the frozen core approximation; rather, an all-electron basis set is used. Scalar relativistic effects are considered using the zeroth-order regular approximation (ZORA). In EDA, bond formation between the interacting fragments is divided into four steps, which can be interpreted in a plausible way. In the first step the fragments, which are calculated with the frozen geometry of the entire molecule, are superimposed without electronic relaxation, yielding the quasiclassical electrostatic attraction ΔEelstat. In the second step the product wave function becomes antisymmetrized and renormalized, which gives the repulsive term ΔEPauli, termed Pauli repulsion. In the third step the molecular orbitals relax to their final form to yield the stabilizing orbital interaction ΔEorb. The latter term can be divided into contributions of orbitals having different symmetries. As it has already been mentioned, the dispersion corrected revPBE-D3 functional is employed, hence the dispersion correction term, ΔEdisp, will be added to the interaction energy (ΔEint) values to describe the total bond energy as
ΔEint = ΔEPauli + ΔEelstat + ΔEorb + ΔEdisp | (1) |
−De = ΔEprep + ΔEint | (2) |
Additionally, the interactions are analyzed using the NCI index.10 In NCI analysis, the mapping of localized binding interactions is done by employing two scalar fields, the electron density (ρ), and the reduced density gradient (s). These two quantities are connected by
(3) |
(4) |
(5) |
To perform such integrations, it is necessary to establish a unique definition of the NCI region. Because the difference between the interacting and noninteracting monomers is directly reflected in the s(ρ) diagram, it is possible to define the NCI region as the points in 3D space with (ρ, s) values lying in the s(ρ) peak. To identify this region, both the monomer and dimer densities must be computed and compared. The lower edge of the monomer s(ρ) curve is splined and all the points of the dimer s(ρ) plot lying below the splined curve are localized in real space. In practice, these integrations are performed numerically, by summation over a cubic grid with 0.1 a.u. increments and cutoffs of ρ = 0.2 a.u. and s = 2.0.
The inclusion of dispersion via Grimme's corrections induces the collapse of some of these isomers, resulting in only four dimers for each metallocene dimer (Fig. 2). Remarkably, regardless of the conformation of the monomers in the initial guess (eclipsed or staggered), dimers have only eclipsed units. When dispersion is involved, structure 1, in which two metallocenes are oriented perpendicular to each other, becomes the most stable one in all cases. Structure 2 is the one described by Bogdanović and Novaković, and is the second most stable form.6 In the higher energy forms, one cyclopentadienyl ring of the first unit interacts perpendicularly (3) or in a parallel manner (4) to the second unit. In the case of ferrocene, while 1 is only 1.7 kcal mol−1 lower in energy than 2, isomers 3 and 4 are 3.0 and 4.7 kcal mol−1 above the global minimum, respectively. The same energy order is noticed for the ruthenocene and osmocene dimers, but the energy differences increase considerably (see Fig. 2).
Dispersion has strong effects on the geometrical parameters. The M⋯M distances without dispersion in isomer 1 are 5.80 (Fe), 5.56 (Ru), and 5.48 (Os) Å. These distances are reduced to 4.84 (Fe), 4.72 (Ru), and 4.64 (Os) Å by applying the D2 correction to PBE, i.e., the intermolecular dispersion reduces the M⋯M distances by almost 1 Å. Note that a reduction in the M⋯M distance from Fe to Os is independent of the dispersion inclusion, indicating a stronger interaction for the osmocene dimer. A similar shortening of the M⋯M distance is perceived in the other three isomers as a consequence of the dispersion inclusion. At the PBE-D2/def2-TZVP level, with inclusion of the ZPE correction, the computed bond dissociation energies (BDE) for 1 without dispersion are 0.7 (Fe), 1.3 (Ru), and 1.0 (Os) kcal mol−1. These values are remarkably increased to 7.5 (Fe), 10.0 (Ru), and 13.7 (Os) kcal mol−1, when dispersion is included.
Temperature also significantly affects the dissociation energies. Fig. 3b indicates that while the ferrocene dimer is stable until 220 K, the ruthenocene and osmocene dimers are stable complexes until 270 and 350 K, respectively. The main reason is that at high temperatures, the contact area between the two units decreases and less compact clusters are obtained. In other words, the ideal gas behavior (no intermolecular interactions, thus no cluster formation) is recovered at high temperatures.
M | Cluster | ΔEelstat | ΔEPauli | ΔEorb | ΔEdisp | ΔEtotal |
---|---|---|---|---|---|---|
The percentage values within the parentheses show the contribution towards the total attractive interactions (ΔEelstat + ΔEorb + ΔEdisp). | ||||||
Fe | 1 | −9.3 (30.0) | 23.3 | −5.3 (17.1) | −16.4 (52.9) | −7.7 |
2 | −6.2 (28.4) | 16.1 | −3.7 (17.0) | −11.9 (54.6) | −5.7 | |
Ru | 1 | −11.1 (31.4) | 26.5 | −6.6 (18.6) | −17.7 (50.0) | −9.0 |
2 | −8.2 (32.2) | 18.8 | −4.6 (18.0) | −12.7 (49.8) | −6.7 | |
Os | 1 | −13.4 (32.1) | 30.2 | −7.5 (17.9) | −20.9 (50.0) | −11.6 |
2 | −9.7 (32.7) | 21.5 | −5.2 (17.5) | −14.8 (49.8) | −8.2 | |
(H2O)2 | −8.8 (59.9) | 9.9 | −4.7 (32.0) | −1.2 (8.2) | −4.8 | |
(C6H6)2 | −2.5 (24.5) | 7.3 | −1.1 (10.8) | −6.6 (64.7) | −2.8 |
Bogdanović and Novaković found an electrostatic complementarity between the two ferrocene units in structure 2.6 This complementarity occurs in a very large area, including the four “puzzle-like” regions of mutual compatibility and recognition. Fig. 4 shows molecular electrostatic potential isosurfaces for the four local minima of ferrocene. It is apparent that this complementarity concept is also applicable to structure 1. However, this does not apply to 3 and 4, where the overlap of the negative regions interrupts the complementarity.
Fig. 4 Molecular electrostatic potential maps of the ferrocene dimers. Isosurfaces are plotted with an isovalue of 0.01. Positive regions are in red and negative regions in blue. |
Charge density difference plots give a precise representation of the electron redistribution upon dimer formation (Fig. 5). Clearly, the electron redistribution in the three dimers is negligible. The most relevant changes are perceived in the density of the C–H bond involved in the contact. While the hydrogen atoms lose density, carbon atoms gain it. So, a very poor polarization contribution is expected (vide infra). Note that dispersion forces always will be attractive even without any charge transfer.
Fig. 5 Electron density differences in isosurfaces (Δρ = 0.0005 a.u.). Negative and positive regions are in red and blue, respectively. |
The strength of the overall contributions of dispersion and electrostatic terms to the interactions can be quantitatively estimated by EDA. Table 1 gives the EDA values at the revPBE-D3/TZ2P//PBE-D2/def2-TZVP level calculated for the two most stable isomers. It becomes apparent that the electrostatic contribution is not the sole term stabilizing the dimers. It contributes ca. 28–32% towards the total attraction. ΔEdisp is found to be the major contributor towards the total attraction with ≈50%. The contribution from ΔEorb is the least towards the total attraction (less than 20%) as a consequence of negligible orbital overlapping between two monomers. Note that the interaction energy gradually increases on moving from Fe to Os. It is also important to mention that monomers do not suffer any significant structural variation, thus preparation energy is negligible. Clearly, EDA shows a delicate balance between dispersion and electrostatic contributions in order to stabilize a metallocene dimer.
Note that although the magnitude of the interaction energy between two metallocene units is quite close to typical H-bond energy, the nature of the bonding between them is quite different. As a reference point, the water dimer is considered (see Table 1). In this H-bonded dimer, the ΔEelstat contribution is found to be the most significant (60%), whereas the contributions from ΔEorb and ΔEdisp terms to the total attractive interaction are 32 and 8.2%, respectively. Therefore, according to EDA, in the water dimer, as in most cases, the H-bonding is mainly electrostatic in nature.28 We have also compared the title cases with the parallel-displaced benzene dimer with C2h symmetry. EDA indicates that the benzene dimer is stabilized mainly by dispersion (65%). So, the nature of bonding of metallocene dimers is similar to that of the benzene dimer, but the interaction energy is higher in the first cases.
As the name suggests, the NCI method of Johnson and co-workers has been specifically developed to reveal non-covalent interactions.10Fig. 6 depicts NCI isosurfaces for isomers 1 and 2 to illustrate the nature of the intermolecular interactions. As is mentioned above, a continuous color-coding scheme based on the second derivative of the electron density is used, where strong attractive interactions are represented in blue, weak attractive interactions in green, and strong repulsive interactions in red. The images of both isomers correspond to a typical dimer stabilized mainly by dispersion where attractive surfaces cover a very large area between both units. So, NCI analysis supports the fact that there is a collective action of small interaction energies of dispersive nature distributed in an ample area between both monomers. Furthermore, the area of the green surface in 1 is larger than that in 2. Integration of the volume of this surface for Fe (see ESI†) gives 87.47 a.u. for the parallel conformer and 104.03 a.u. for the perpendicular one. The charges involved within these surfaces are also larger for the orthogonal conformation (1.78 |e| in the parallel case vs. 1.90 |e| for the orthogonal one). These values are in agreement with the relative weight of dispersive and electrostatic contributions in the EDA partition.
Footnote |
† Electronic supplementary information (ESI) available: Optimized Cartesian coordinates for all species. See DOI: 10.1039/c5cp05956a |
This journal is © the Owner Societies 2016 |