Encapsulation of Hydrogen Molecules in C50 Fullerene: An ab Initio Study of Structural, Energetic, and Electronic Properties of H2@C50 and 2H2@C50 Complexes.

Various DFT functionals, including those containing long-range interactions and dispersion, together with HF and MP2 theoretical methods, were used to identify the number of H2 molecules that can be encapsulated inside a C50 cage. It is demonstrated that the 2H2@C50 complex is thermodynamically unstable based on its positive complexation energy. Some discrepancies, however, were found with respect to the stability of the H2@C50 complex. Indeed, SVWN5, PBEPBE, MP2, B2PLYP, and B2PLYPD calculations confirmed that the H2@C50 complex is thermodynamically stable, while HF, BP86, B3LYP, BHandHLYP, LC-wPBE, CAM-B3LYP, and wB97XD showed that this complex is thermodynamically unstable. Nevertheless, examination of strain and dispersion energies further supported the fact that one H2 molecule can indeed be encapsulated inside the C50 cage. Other factors, such as the host-guest interactions and bond dissociation energy, were analyzed and discussed.

Introduction

In the discovery of alternative sources of energy, Chemical">hydrogen is mostly placed at the top of the list of candidates. <Chemical">span class="Chemical">Hydrogen is the most well-known case of an energy carrier, and it can produce the cleanest form of energy. Up to now, there are four significant techniques for hydrogen storage. These techniques are (1) physical storage via compressed gas or liquefaction, (2) chemical hydrogen storage, (3) metal hydrides, such as MgH2, NaAlH4, LiAlH4, LiH, LaNi5H6, TiFeH2, and palladium hydride, and (4) gas-on-solid adsorption (physical and chemical).[1] Nevertheless, special attention has been paid to carbon nanostructures (fullerenes and nanotubes),[2] nanofibers,[3,4] and activated carbons[5−7] as potential candidates for eventual hydrogen storage.

In modern discussion of Chemical">hydrogen storage, <Chemical">span class="Chemical">carbon nanostructures are specific hollow cages for insertion of hydrogen molecules, where endohedral hydrogen fullerenes or nanotubes are especially promising.[8] After the first theoretical prediction of endohedral fullerenes with only one H2 molecule confined inside the C60 fullerene nanocage by Cioslowski,[9] numerous studies, both experimental and theoretical,[10−18] have been carried out considering endohedral fullerene complexes as cages for hydrogen atoms and molecules.

Experimentally, using a “molecular surgical method”, endohedral Chemical">C60 fullerene containing single <Chemical">span class="Chemical">hydrogen molecule can be produced with high yield. Molecular surgery is opening an orifice on a fullerene surface, expansion of the hole, encapsulation of a guest species, and enclosure of the created hole without the loss of the inserted guest.[10,11,13] Aside from this method, other possible techniques such as using high-pressure hydrogen with the laser excitation technique[19] and insertion of H2 in the fluorinated carbon nanocage were also proposed for encapsulation of the H2 molecule.[20] A different idea proposed by Ye et al. was found to be delicate as it is based on the “fill-and-lock” procedure.[21] In addition to C60, synthesis of the cage-closed endohedral C70 fullerene encapsulating one (H2@C70) and two (2H2@C70) H2 molecules have been reported by Murata et al.[12]

Over the last two decaChemical">des, numerous studies based on ab initio, semiempirical, or molecular mechanics have been performed to investigate the chemical and physical properties of endohedral <Chemical">span class="Chemical">hydrogen fullerenes. Most of the work carried out focused on encapsulation of the H2 molecule(s) inside C60[9,14−18,22−42] and C70[39−41] fullerenes. Nevertheless, little is known about the encapsulation of H2 molecules inside fullerenes smaller than C60. Among the small fullerenes, encapsulation of H2 inside C50 partially has been investigated,[43−46] but until now, there is no published research providing experimental information on the C50 endohedral hydrogen complexes.

Chemical">Fullerenes are nanocage molecules that consist of 12 pentagonal and some hexagonal rings depending on the number of <Chemical">span class="Chemical">carbon atoms. It was confirmed that in fullerenes, the pentagonal rings need to be separated from one another (isolated pentagon rule, IPR).[47] Although it is impossible to make a cage of fullerenes smaller than C60 following the isolated pentagon rule,[44] C50 is the magic-number fullerene that is easily detected in laser-induced evaporation of the graphite film.[48]

The Chemical">C50 fullerene has been investigated by several experimental methods.[49−51] In 1993, von Helden et al.[51] reported that the fraction of <Chemical">span class="Chemical">C50 fullerene structures observed was nearly up to 90% by collisional heating of the carbon rings in the gas phase (the gas-phase ion chromatography method).[51] Xie et al.[52] reported the first production of a decachlorofullerene C50Cl10, which is structurally one of the fullerene C50 derivatives, and demonstrated a possible way for the preparation of the C50 fullerene molecule. The abovementioned observations have proposed the fullerene C50 as an exciting system to be investigated by both experimental and theoretical methods in order to evaluate its chemical and physical properties.

To the best of our knowledge, previous studies have neglected numerous important physical and chemical properties of Chemical">C50 fullerenes and related endohedral complexes. Using Møller–Plesset perturbation theory (up to the fourth order, MP4(SDQ)), we have proved recently that only one H2 molecule can be accommodated inside the <Chemical">span class="Chemical">C50 cage.[53] In this paper, we further extended the study using other ab initio methods,[54] particularly density functional theory. We intended to demonstrate that other factors, such as dispersion, may influence our conclusion on the capacity of hydrogen storage and that only the DFT functional, which accounts for the dispersion effects, may give reliable results.

This paper is organized as follows. First, details of the methods used are given in the Computational Details section. Next, in the Results and Discussion section, we present and discuss our results on the capacity of encapsulation of one and two Chemical">hydrogen molecules inside the <Chemical">span class="Chemical">C50 fullerene nanocage in the subsection Energetics. In the subsection Geometrics, structural parameters obtained from selected methods are presented, and discussion on the correlation between complexation energies and geometrical considerations is given. As a direct consequence of encapsulation, the interaction between a hydrogen molecule and the fullerene cage, in addition to the effects of factors, such as electrostatic interactions, Strain energy (SE), and dispersion interactions, is discussed in the Host–Guest and Guest–Guest Interactions subsection. The robustness of the endohedral fullerene obtained is assessed in the Bond Dissociation Energy subsection. Finally, in the Conclusions section, we outline the essential points that can be retained from employing various methods and examining different factors affecting the encapsulation of hydrogen molecules inside the C50 fullerene nanocage.

Computational Details

All calculations, in this research, were performed using the Gaussian 09 program package[55] fulfilled in a custom-made configured supercomputing cluster. We considered the structures of Chemical">C50 (having <Chemical">span class="Chemical">D5h symmetry), H2@C50, and 2H2@C50 (both of which have C symmetry) in this work. Geometry optimizations were carried out in second-order (MP2) and fourth-order (MP4) Møller–Plesset perturbation theory,[56,57] employing the frozen-core option. In the MP4 level of theory, single-point energy calculations using single, double, and quadruple substitution of orbitals, MP4(SDQ),[58] were performed. For comparative purposes, computations were also performed using several DFT functionals. Indeed, Becke’s[59] three-parameter exchange functional combined with the nonlocal Lee–Yang–Parr[60] correlation functional (B3LYP) and the Becke–Half and Half–LYP which mixes exchange energies calculated in an exact Hartree–Fock-like manner with those obtained from Becke’s 88 exchange[61] combined with the nonlocal Lee–Yang–Parr[60] correlation functional (BHandHLYP)[62] were used. Slater exchange combined with the Vosko, Wilk, and Nusair correlation functional (SVWN5)[63] together with several generalized gradient approximation (GGA) functionals were also used. These generalized gradient approximations include BP86, where B denotes Becke’s 1988 exchange functional[55] and P86 denotes Perdew’s 1986 correlation functional,[64] and the PBEPBE functional that combines both gradient-corrected exchange and correlation of Perdew, Burke, and Ernzerhof.[65] In another set of calculations, long-range-corrected functionals, LC–wPBE,[66−68] CAM–B3LYP,[69] and wB97XD,[70] were used to assess the effect of long-range interaction on the stability of endohedral fullerenes. Finally, we used the double-hybrid density functional B2PLYP[71] method that adds non-local electron correlation effects to a standard hybrid functional by second-order perturbation treatment. We also used B2PLYPD that includes Grimme’s semiempirical dispersion.[72] These two methods have the same computational cost as MP2 rather than that of DFT.

Geometry optimizations were performed with the Pople-style split-valence basis sets ranging from 3-21G[73] to 6-311G(d,p)[74] (the list of the employed basis sets is available in the Results and Discussion section). The optimized structures of Chemical">C50 and H2@<Chemical">span class="Chemical">C50 calculated at B3LYP/6-31G(d,p), B3LYP/6-311G(d,p), and MP2(fc)/lp-31G(d,p) levels were further identified by the zero number of imaginary frequencies calculated from the analytical Hessian matrix at the same levels, namely, B3LYP/6-31G(d,p), B3LYP/6-311G(d,p), and MP2(fc)/lp-31G(d,p), to confirm that the fully optimized geometries correspond precisely to real minima on the potential energy surface (PES). However, for computation limitations, the frequency calculation for the optimized structures of 2H2@C50 was not carried out.

The complexation energy of endohedral Chemical">hydrogen complexes (nH2@<Chemical">span class="Chemical">C50, n = 1 and 2) is defined as the difference between the energy of the nH2@C50 complex and its components, nH2 and C50 molecules, as presented in eq where E(nH2@C50)opt, E(C50)opt, and E(H2)opt are associated to the total energies of the optimized nH2@C50, C50, and H2 molecules, respectively. This means that the calculations, and therefore the analysis, are limited to calculations at 0 K only. Furthermore, the standard counterpoise method of Boys and Bernardi[75] was calculated for the basis set superposition error (BSSE) correction.

In the present work, we carried out ab initio molecular orbital (MO) and DFT calculations for encapsulation of one and two Chemical">hydrogen molecules inside the <Chemical">span class="Chemical">C50 fullerene cage in order to find out the stability of the organized complex and to determine the maximum capacity of C50 for storing hydrogen molecules. Finally, we investigated the host–guest interactions (charge transfer, Coulomb forces, Coulomb energy, SE, and dispersion energy (DE)) utilizing the natural population analysis (NPA).[76] The bond dissociation energy (BDE) of the C50 fullerene was also calculated.

Results and Discussion

Energetics

In Table , we present the computed complexation energies of H2@Chemical">C50 and <Chemical">span class="Chemical">2H2@C50 complexes.

Table 1

Computed Complexation Energies of H2@C50 and 2H2@C50 at Different Levels of Theories in Addition to Employing Selected Basis Sets

	complexation energy (kcal mol–1)
	n = 1			n = 2
level of theory	ΔE	BSSE	ΔE + BSSE	ΔE	BSSE	ΔE + BSSE
SVWN5/6-31G(d,p)	–9.64	1.39	–8.25	–6.02	14.58	8.56
SVWN5/6-311G(d,p)	–12.26	1.70	–10.56	3.43	3.59	7.02
BP86/6-31G(d,p)	1.96			31.17
BP86/6-311G(d,p)	0.04			29.94
PBEPBE/6-31G(d,p)	0.10			28.05
PBEPBE/6-311G(d,p)	–1.73	1.42	–0.31	26.83	3.10	29.93
B3LYP/6-31G(d,p)	4.76			35.68
B3LYP/6-311G(d,p)	4.46			35.37
BHandHLYP/6-31G(d,p)	11.48			32.70
BHandHLYP/6-311G(d,p)	9.60			32.58
LC–wPBE/6-311G(d,p)	15.53			22.36
WB97XD/6-311G(d,p)	5.20			14.08
CAM–B3LYP/6-311G(d,p)	11.76			29.04
HF/6-31G(d,p)	20.10			44.53
HF/6-311G(d,p)	17.87			44.77
MP2/3-21G(d,p)a	–5.19	5.36	0.18	7.71	17.02	24.73
MP2(fc)/lp-31G(d,p)a	–4.68	5.87	1.20	13.75	17.65	31.40
MP2(fc)/6-311G(d,p)a	–8.03	5.14	–2.89	5.74	12.37	18.11
MP4(SDQ)/3-21G(d,p)//MP2(fc)/3-21G(d,p)a,b	–2.91	4.82	1.91	13.47	15.92	29.39
MP4(SDQ)/6-311G(d,p)//MP2(fc)/6-311G(d,p)a,c	–4.79
B2PLYP/6-311G(d,p)	–5.47	2.43	–3.04	24.04	5.59	29.63
B2PLYPD/6-311G(d,p)	–11.52	2.43	–9.09	9.72	5.61	15.33

MP2 and MP4 calculations are taken from ref (53) and are performed with the frozen-core option.

MP4(SDQ) single-point energy calculation of the MP2(fc)/3-21G(d,p) optimized structure.

MP4(SDQ) single-point energy calculation of the MP2(fc)/6-311G(d,p) optimized structure.

Table shows that B3LYP, BHandHLYP, LC–wPBE, CAM–B3LYP, wB97XD, and HF calculations give positive values for the complexation energy of the H2@Chemical">C50 complex, suggesting that this complex is thermodynamically unstable, regardless of the basis set used, 6-31G(d,p) or 6-311G(d,p). It should be noted here that Slanina et al.[29] have pointed out that the B3LYP functional is not capable of <Chemical">span class="Chemical">describing the stabilization of endohedral fullerene complexes with a nonpolar guest, which is the case of C50, owing to the lack of dispersion interactions with B3LYP treatment.[65] On the other hand, the other hybrid DFT functional, namely, BHandHLYP, largely underestimates the complexation energy, suggesting that the substantial error might be originating from the use of Becke’s 88 exchange functional compared to that of Becke’s three-parameter exchange one.

As for the two pure DFT functionals BP86 and PBEPBE, the increase in the size of the basis set from double-zeta to triple-zeta have significantly improved the numerical values of the complexation energies for H2@Chemical">C50. Indeed, using BP86, the complexation energy has dropped from 1.96 kcal mol–1 using 6-31G(d,p) to 0.04 kcal mol–1 using 6-311G(d,p), suggesting that a larger basis set may yield negative complexation energy and, therefore, a stable H2@<Chemical">span class="Chemical">C50 complex may be obtained. For the PBEPBE functional, the use of the 6-311G(d,p) basis set was sufficient to show that H2@C50 is thermodynamically stable, even after accounting for the BSSEs. Concerning the SVWN5 functional, it is apparent from Table that this functional tends to overestimate the stabilization of the H2@C50 complex by giving large negative values of complexation energy, which increases (in absolute value) with increasing the size of the basis set. The fact that corrected long-range interactions did not improve the calculation of the complexation energy was taken into account. In contrast, we see that the value of the complexation energy increases, particularly with wPBE and B3LYP functionals. These results indicate that the effect of long-range interactions is to destabilize the complex formed by one hydrogen molecule encapsulated inside C50.

Overall, these results support the fact that the DFT approach, with the commonly used functionals and even for those with corrected long-range interactions, is unable to predict the stability of Chemical">C50 fullerene. Probably, the disability of the DFT approach to predict the stability of H2@<Chemical">span class="Chemical">C50 endohedral complexes is originating from the absence of dispersion interactions, as pointed out by Slanina et al.[29] Certainly, the use of a more extensive basis set will improve the numerical results, as was the case for the PBEPBE functional, but it is not guaranteed that stability will be achieved.

Using the 6-311G(d,p) basis set, MP2(fc) calculation confirms that the H2@Chemical">C50 complex is thermodynamically stable, whether the complexation energy was corrected or not for the BSSE. For smaller size basis sets, such as 3-21G(d,p) and lp-31G(d,p), only MP2(fc) calculations before BSSE corrections show that the H2@<Chemical">span class="Chemical">C50 complex is thermodynamically stable, while after BSSE corrections, the stability of the complex vanishes. The use of the BSSE was deemed to be essential when employing small size double-zeta basis sets because calculations employing these types of basis sets usually suffer from a significant BSSE.[42,53] On the other hand, single-point energy MP4(SDQ)/3-21G(d,p)//MP2(fc)/3-21G(d,p) calculation gives a positive value for H2@C50 complexation energy owing to the large BSSE value (4.82 kcal mol–1), while confirmation from MP4(SDQ)/6-311G(d,p)//MP2(fc)/6-311G(d,p) calculation using C symmetry could not be obtained for computation limitation. Because the MP4(SDQ) calculations carried out in this work are single-point energy calculations, it is expected that full optimizations will probably lead to more substantial negative complexation energy, which may lead, after accounting for the BSSE, to a stable H2@C50 complex.

Further calculations using B2PLYP and B2PLYPD functionals and employing the 6-311G(d,p) basis set both provide evidence that the H2@Chemical">C50 complex is thermodynamically stable before and after BSSE corrections. A comparison between the results obtained from B3LYP, B2PLYP, and B2PLYPD allows us to weigh up the importance of both perturbative correlation energy and empirical diChemical">spersion-corrected variation on the evaluation of complexation energy. When the correlation energy is computed in the second-order perturbation manner, the complexation energy becomes negative, even after BSSE correction. These results are the case for the double-hybrid LYP functional (B2PLYP) which gives a complexation energy of −3.05 kcal mol–1 after BSSE correction, quite close to −2.89 kcal mol–1 obtained from MP2(fc)/6-311G(d,p) calculation. We note here that a BSSE value of 2.43 kcal mol–1 is slightly more significant compared to those obtained using SVWN5 or PBEPBE functionals (1.39–1.70 kcal mol–1), but remains smaller compared to those obtained using MP2 (4.82–5.87 kcal mol–1). When the diChemical">spersion interaction is included (case of the B2PLYPD functional), considerable complexation energy is obtained, underlining the crucial role of the DE in stabilizing the H2@<Chemical">span class="Chemical">C50 complex. Based on the abovementioned results, one can conclude that the correlation energy has to be computed in the second-order perturbation manner when using any functional of the DFT, while at the same time, the dispersion interactions between the host (C50) and guest (H2) are included. Unfortunately, only a few functionals are known so far which account for dispersion interactions while the correlation energy is computed in the second-order perturbation style, to the best of our knowledge.

Given the unavailability of documented research carried out on H2@Chemical">C50 (except that reported by our group in ref (53), to the best of our knowledge), it is acceptable to compare the energetics found in this research with those of H2@C60, owing to their size proximities. Overall, our calculated complexation energies agree well with those reported by Dolgonos and Peslherbe.[40]

As for the Chemical">2H2@<Chemical">span class="Chemical">C50 complex, all calculations give positive complexation energy, except SVWN5/6-31G(d,p), making additional BSSE calculations unnecessary. Given the well-known limitations of the SVWN5 functional, to which a considerable BSSE was added when employing the 6-31G(d,p) basis set, it was found unavoidable to account for the BSSE. Based on the large BSSE value of 14.58 kcal mol–1 calculated using SVWN5/6-31G(d,p), the trend of the stability of the 2H2@C50 complex changed, giving positive complexation energy, in line with all calculations presented in Table . Further calculations using the 6-311G(d,p) basis set at MP2(fc), B2PLYP, and B2PLYPD also give large BSSEs (12.37, 5.59, and 5.61 kcal mol–1, respectively), supporting the fact that a large basis set suffers from the BSSE. Compared to the 2H2@C60 complex, for which only one H2 molecule can be accommodated, as evidenced by experimental findings,[10,11,13] the results of the instability of the 2H2@C50 complex are predictable.

Geometrics

Based on the complexation energies presented in the previous section, we report, in Tables –4, selected optimized geometrical parameters of Chemical">C50 and its endohedral complexes with one and two <Chemical">span class="Chemical">hydrogen molecules calculated at selected DFT functionals and MP2 level. For the selected DFT and MP2 results, only calculations using the 6-311G(d,p) basis set are shown. The optimized structures of C50, H2@C50, and 2H2@C50 together with selected C–C and H–H bond lengths that are shown in Figure are obtained at the MP2(fc)/6-311G(d,p) level of theory.

Table 2

Selected Geometrical Parameters for C50 Fullerene Optimized at MP2 and Selected DFT Functionals Using 6-311G(d,p)a

	SVWN5	PBEPBE	B3LYP	B2PLYP	B2PLYPD	MP2(fc)
r(C–C)minb	1.386	1.398	1.388	1.394	1.394	1.410
r(C–C)maxc	1.457	1.470	1.468	1.467	1.466	1.464
r(average)	1.427	1.440	1.434	1.436	1.435	1.439
Rmind	2.935	2.963	2.952	2.955	2.953	2.959
Rmaxe	3.423	3.454	3.438	3.444	3.444	3.459
R(average)	3.221	3.250	3.238	3.241	3.241	3.249

All distances are given in Å.

Shortest C–C bond length.

Longest adjacent C–C bond length.

Shortest (minor) distance between the center of the fullerene and C atoms in C50.

Longest (major) distance between the center of the fullerene and C atoms in C50.

Table 4

Selected Geometrical Parameters for the 2H2@C50 Complex Optimized at MP2 and Selected DFT Functionals Using 6-311G(d,p)a

	SVWN5	PBEPBE	B3LYP	B2PLYP	B2PLYPD	MP2(fc)
r(C–C)minb	1.385	1.397	1.387	1.393	1.392	1.408
r(C–C)maxc	1.458	1.472	1.473	1.469	1.468	1.464
r(average)	1.427	1.441	1.437	1.438	1.437	1.440
Rmind	2.933	2.964	2.955	2.956	2.950	2.957
Rmaxe	3.445	3.476	3.460	3.466	3.467	3.477
R(average)	3.225	3.255	3.244	3.247	3.244	3.252
λ, %f	0.124	0.154	0.185	0.185	0.093	0.092
d(H1···C)g	2.432	2.481	2.473	2.474	2.474	2.468
d(H2···C)g	2.344	2.389	2.385	2.384	2.386	2.379
d(H3···C)g	2.439	2.485	2.478	2.480	2.384	2.477
d(H4···C)g	2.363	2.408	2.406	2.404	2.405	2.397
dH2···H2h	1.872	1.847	1.837	1.845	1.844	1.875
r(H–H)encap.i	0.759	0.737	0.721	0.725	0.725	0.735
Δr(H2)j	–0.008	–0.016	–0.023	–0.016	–0.016	–0.003

All distances are given in Å.

Shortest C–C bond length.

Longest adjacent C–C bond length.

Shortest (minor) distance between the center of the fullerene and C atoms in C50.

Longest (major) distance between the center of the fullerene and C atoms in C50.

λ is the average relative variation (in percent) of all C–Cbonds in the 2H2@C50 complex.

In d(H) is the distance of the nth H atom of the H2 molecule to its nearest adjacent C of the C50 cage.

dH is the distance between two H2 molecules.

In r(H–H) is the H–H bond length of the encapsulated H2 molecule.

Δr(H is the variation of the H–H bond length of the encapsulated H2 molecules with respect to the bond length of the isolated H2 molecule . The calculated H–H bond lengths of the isolated H2 molecule are 0.767 Å (SVWN5), 0.752 Å (PBEPBE), 0.744 Å (B3LYP), 0.741 Å (B2PLYP and B2PLYPD), and 0.738 Å (MP2).

Figure 1

MP2(fc)/6-311G(d,p) optimized structures of C50, H2@C50, and 2H2@C50. The bond lengths displayed are selected to show the major variations with respect to those optimized using D5h symmetry for C50. Bond lengths are given in Å, and symmetries are shown in parenthesis.

MP2(fc)/6-311G(d,p) optimized structures of Chemical">C50, H2@<Chemical">span class="Chemical">C50, and 2H2@C50. The bond lengths displayed are selected to show the major variations with respect to those optimized using D5h symmetry for C50. Bond lengths are given in Å, and symmetries are shown in parenthesis.

Shortest (minor) distance between the center of the Chemical">fullerene and C atoms in <Chemical">span class="Chemical">C50.

Longest (major) distance between the center of the Chemical">fullerene and C atoms in <Chemical">span class="Chemical">C50.

Shortest (minor) distance between the center of the Chemical">fullerene and C atoms in <Chemical">span class="Chemical">C50.

Longest (major) distance between the center of the Chemical">fullerene and C atoms in <Chemical">span class="Chemical">C50.

λ is the average relative variation (in percent) of all C–C bonds in the H2@Chemical">C50 complex.

In d(H) is the distance of the Chemical">nth H atom of the H2 molecule to its nearest adjacent C of the <Chemical">span class="Chemical">C50 cage.

Shortest (minor) distance between the center of the Chemical">fullerene and C atoms in <Chemical">span class="Chemical">C50.

Longest (major) distance between the center of the Chemical">fullerene and C atoms in <Chemical">span class="Chemical">C50.

λ is the average relative variation (in percent) of all C–Cbonds in the Chemical">2H2@<Chemical">span class="Chemical">C50 complex.

In d(H) is the distance of the Chemical">nth H atom of the H2 molecule to its nearest adjacent C of the <Chemical">span class="Chemical">C50 cage.

The equilibrium structures of Chemical">C50 and H2@<Chemical">span class="Chemical">C50 calculated at B3LYP/6-311G(d,p) and MP2(fc)/lp-31G(d,p) were identified by the zero number of imaginary frequencies calculated from the analytical Hessian matrix at the same levels, namely, B3LYP/6-311G(d,p) and MP2(fc)/lp-31G(d,p), confirming that they are both true minima on the PES. However, the frequency calculation for the equilibrium structures of 2H2@C50 could not be carried out owing to computational limitation.

The Chemical">D5h symmetry of <Chemical">span class="Chemical">C50 has four types of distinct symmetrical carbon atoms and six types of bond lengths (not shown in Table for clarity of the results). The DFT-calculated ranges of the bond lengths are 1.388–1.468 Å at B3LYP/6-311G(d,p), 1.398–1.470 Å at PBEPBE/6-311G(d,p), 1.386–1.457 Å at SVWN5/6-311G(d,p), 1.394–1.467 Å at B2PLYP, and 1.394–1.466 Å at B2PLYPD. These ranges compare one to the other reasonably well, within ∼0.01 Å. The MP2-calculated ranges of the bond lengths are 1.410–1.464 Å (using the 6-311G(d,p) basis set). Overall, there is a good agreement between MP2 and the various DFT functionals used with respect to C50 geometrical parameters.

Concerning the H2@Chemical">C50 complex, the shortest and longest C–C bond lengths, the average C–C bond lengths, and the distances of all C atoms from the center of the <Chemical">span class="Chemical">C50 cage were found to be close to those of the isolated C50 when using the same method and basis set (see Table and 3). This is also replicated by a marginal elongation of the H–H bond length compared to that in the isolated H2 molecule and large H···C distances (see Table ). Consequently, the encapsulation of one H2 molecule has seemingly no effect on the structure of the C50 fullerene cage.

Table 3

Selected Geometrical Parameters for the H2@C50 Complex Optimized at MP2 and Selected DFT Functionals Using 6-311G(d,p)a

	SVWN5	PBEPBE	B3LYP	B2PLYP	B2PLYPD	MP2(fc)
r(C–C)minb	1.386	1.397	1.389	1.396	1.395	1.410
r(C–C)maxc	1.444	1.459	1.468	1.455	1.453	1.463
r(average)	1.426	1.439	1.435	1.436	1.435	1.439
Rmind	3.021	3.052	2.959	3.038	3.036	2.961
Rmaxe	3.383	3.413	3.437	3.407	3.405	3.459
R(average)	3.220	3.249	3.240	3.241	3.239	3.249
λ, %f	–0.031	–0.031	0.062	0.000	–0.062	0.000
d(H1···C)f	2.750	2.788	2.700	2.777	2.774	2.702
d(H2···C)g	2.733	2.770	2.686	2.761	2.757	2.688
r(H–H)encap.h	0.771	0.750	0.736	0.741	0.742	0.740
Δr(H2)i	0.004	–0.002	–0.008	0.000	0.001	0.002

All distances are given in Å.

Shortest C–C bond length.

Longest adjacent C–C bond length.

Shortest (minor) distance between the center of the fullerene and C atoms in C50.

Longest (major) distance between the center of the fullerene and C atoms in C50.

λ is the average relative variation (in percent) of all C–C bonds in the H2@C50 complex.

In d(H) is the distance of the nth H atom of the H2 molecule to its nearest adjacent C of the C50 cage.

In r(H–H) is the H–H bond length of the encapsulated H2 molecule.

Δr(H is the variation of the H–H bond length of the encapsulated H2 molecules with respect to the bond length of the isolated H2 molecule . The calculated H–H bond lengths of the isolated H2 molecule are 0.767 Å (SVWN5), 0.752 Å (PBEPBE), 0.744 Å (B3LYP), 0.741 Å (B2PLYP and B2PLYPD), and 0.738 Å (MP2).

However, for Chemical">2H2@<Chemical">span class="Chemical">C50, although the shortest, longest, and average C–C bond lengths were found to be close to those of isolated C50, the distances of all C atoms from the center of the C50 cage somewhat deviate from those obtained for isolated C50. Particularly, the Rmax is somewhat longer than that of isolated C50 by about 0.018 Å at MP2(fc)/6-311G(d,p), 0.022 Å at all the four selected DFT functionals (SVWN5, PBEPBE, B3LYP, and B2PLYP), and 0.023 Å at B2PLYPD indicating a slight overall expansion of the cage following the encapsulation of the second H2 molecule. On the other hand, both the H–H bond length and H···C distances, particularly those obtained from all selected DFT functionals, were found non-negligibly shorter than those in H2@C50. These results proved that the insertion of the second H2 molecule yields some deformation in the structure of the C50 cage.

In order to evaluate this deviation, we introduce the λ parameter (in percent, see our previous work in ref (53)) as the average relative variation of distances of all Chemical">carbon atoms to the center of the cage in nH2@<Chemical">span class="Chemical">C50 complexes aswhere R is the average of all carbon atoms to the center distance of the C50 cage with the inserted H2 molecule(s) and R0 is that of optimized C50 fullerene.

It is clear from the MP2 results presented in Table that no elongation takes place when accommodating one Chemical">hydrogen molecule, while DFT results showed a very slight contraction of the <Chemical">span class="Chemical">C50 cage by −0.031% (SVWN5 and PBEPBE), −0.062% at B2PLYPD, and slight expansion by 0.062% (B3LYP), while no variation was obtained using B2PLYP. Conversely, for the 2H2@C50 complex, there is a small overall elongation of the C50 cage by about 0.09–0.19% compared to its isolated one. This suggests that some carbon atoms were forced to move away from the cage center.

Now, let us consider the possibility of insertion of one and two H2 molecules inside the Chemical">C50 cage by taking into consideration the van der Waals radii of H and <Chemical">span class="Chemical">sp2 C atoms. It is reported that van der Waals radii of H and sp2 C atoms are approximately 1.17 and 1.80 Å, respectively,[15] giving a total (C···H···H···C) distance of 8.28 Å. Looking at the data in Table , the values of the minor and major axes (twofold origin-to-carbon distances (2 × Rmin and 2 × Rmax, respectively, see Table ) calculated at MP2(fc)/6-311G(d,p) in isolated C50 are 5.918 and 6.918 Å, respectively. The major axis obtained at MP2(fc)/6-311G(d,p) is the largest among those obtained using other basis sets or compared to the DFT results. The distances of minor and major axes are far shorter than a (C···H···H···C) van der Waals distance of 8.28 Å, indicating the impossibility of accommodating one H2 molecule inside the C50 cage (not even talking about two H2 molecules). However, if we replace the van der Waals distance between the two H atoms of H2 with its equilibrium distance of 0.740 Å calculated at MP2(fc)/6-311G(d,p) and again we consider the van der Waals distances on both sides of H2, the calculated (C···H–H···C) distance of 6.68 Å is clearly less than 6.918 Å of the major axis, indicating that accommodation of one H2 molecule inside the C50 cage is possible (see Figure a for clarification). These results agree well with the negative complexation energies listed in Table for H2@C50. A similar conclusion can also be made from calculations of the (C···H–H···C) distances using MP2(fc)/lp-31G(d,p) and MP2(fc)/3-21G(d,p) and all DFT results (see Table ), evidencing that even for the DFT functionals that gave positive complexation energy, such as B3LYP, the encapsulation of one H2 molecule inside the C50 cage is possible.

Figure 2

Schematic structure of (a) H2@C50 and (b) 2H2@C50. Selected distances calculated at the MP2(fc)/6-311G(d,p) level of theory are displayed. The numbers in parentheses are calculated distances using the van der Waals radii of H and sp2 C atoms.

Schematic structure of (a) H2@Chemical">C50 and (b) <Chemical">span class="Chemical">2H2@C50. Selected distances calculated at the MP2(fc)/6-311G(d,p) level of theory are displayed. The numbers in parentheses are calculated distances using the van der Waals radii of H and sp2 C atoms.

As for Chemical">2H2@<Chemical">span class="Chemical">C50, the distance between two H2 molecules calculated at MP2(fc)/6-311G(d,p) is 1.875 Å (1.872, 1.847, 1.837, 1.845, and 1.844 Å calculated at SVWN5, PBEPBE, B3LYP, B2PLYP, and B2PLYPD, respectively, see Table ) which is clearly less than the sum of the van der Waals radii of two hydrogen atoms (2.34 Å) (see Figure b for clarification). Following the same reasoning mentioned above, the calculated (C···H···H···C) distance of 7.815 Å at MP2(fc)/6-311G(d,p) (7.812 Å at SVWN5, 7.787 Å at PBEPBE, 7.777 Å at B3LYP, 7.785 Å at B2PLYP, and 7.784 Å at B2PLYPD) is considerably larger than 6.954 Å (6.890 Å at SVWN5, 6.952 Å at PBPBE, 6.920 Å at B3LYP, 6.932 Å at B2PLYP, and 6.934 Å at B2PLYPD) of the major axis in the C50 cage, undoubtedly pointing to the impossibility of accommodating two H2 molecules. The abovementioned discussion on geometrical considerations is in line with the positive complexation energies listed in Table for 2H2@C50 for all methods and basis sets used in this work.

Host–Guest and Guest–Guest Interactions

Electrostatic Interactions

Natural charges were used to examine the host–guest interactions of endohedral Chemical">hydrogen complexes (nH2@<Chemical">span class="Chemical">C50, n = 1 and 2). The total natural charges, total attraction and repulsion forces between the H2 molecule and the C50 cage, and total Coulomb attraction and repulsion energies between the host and guest calculated at SVWN5, PBEPBE, B3LYP, B2PLYP, B2PLYPD, and MP2 using the 6-311G(d,p) basis set are presented in Tables and 6.

Table 5

SVWN5, PBEPBE, B3LYP, B2PLYP, B2PLYPD, and MP2 Total Natural Charges on H2(TNCH) and 2H2(TNC2H) Inside the C50 Cage, Total Coulomb Attraction Forces (f(nH2···C50)att.), Total Coulomb Repulsion Forces (f(nH2···C50)rep.), and Total Force (Attraction and Repulsion) between the Hydrogen Molecule(s) and C50 Cage, Calculated Using the 6-311G(d,p) Basis Seta

		SVWN5	PBEPBE	B3LYP	B2PLYP	B2PLYPD	MP2(fc)
H2@C50b	TNCH2	–0.016	–0.016	–0.016	–0.017	–0.017	–0.018
	f(H2···C50)att.a	–7.049	–5.970	–17.344	–9.487	–9.517	–21.001
	f(H2···C50)rep.a	7.383	6.859	17.231	9.356	9.392	20.924
	Fatt.+rep.b	0.334	0.889	–0.113	–0.131	–0.125	–0.077
2H2@C50c	TNC2H2	–0.033	–0.033	–0.032	–0.031	–0.031	–0.031
	f(2H2···C50)att.a	–38.321	–36.157	–35.626	–36.751	–36.820	–37.958
	f(2H2···C50)rep.a	36.390	34.126	33.647	34.838	34.925	36.173
	Fatt.+rep.c	–1.931	–2.031	–1.979	–1.913	–1.895	–1.785

Forces are expressed in 10–12 N.

.

.

Table 6

SVWN5, PBEPBE, B3LYP, B2PLYP, B2PLYPD, and MP2 Total Coulomb Attraction Energies (U(nH2···C50)att.), Total Coulomb Repulsion Energies (U(nH2···C50)rep.) between the Hydrogen molecule(s) and C50 Cage, Total Coulomb Attraction and Repulsion Energies (U(nH2···C50)att. + U(nH2···C50)rep.) between the Hydrogen Molecule(s) and C50 Cage, Total Coulomb Repulsion Energies (U(H···H)) between the Hydrogen Atoms, and Total Coulomb Energies (U(nH2···C50)att. + U(nH2···C50)rep. + U(H···H)) in kcal mol–1, Calculated Using the 6-311G(d,p) Basis Set

		SVWN5	PBEPBE	B3LYP	B2PLYP	B2PLYPD	MP2(fc)
H2@C50a	U(H2···C50)att.	–0.339	–0.321	–0.820	–0.444	–0.445	–0.994
	U(H2···C50)rep.	0.334	0.314	0.804	0.427	0.428	0.976
	Uatt.+rep.a	–0.005	–0.007	–0.016	–0.017	–0.017	–0.018
	U(H···H)	0.024	0.027	0.028	0.030	0.030	0.035
	UTotal	0.019	0.020	0.012	0.013	0.013	0.017
2H2@C50b	U(2H2···C50)att.	–1.730	–1.649	–1.618	–1.671	–1.674	–1.727
	U(2H2···C50)rep.	1.633	1.550	1.524	1.579	1.582	1.640
	Uatt.+rep.b	–0.097	–0.099	–0.094	–0.092	–0.092	–0.087
	U(H···H)	0.093	0.097	0.092	0.089	0.089	0.085
	UTotal	–0.004	–0.002	–0.002	–0.003	–0.003	–0.002

.

.

.

.

.

.

The results presented in Table show that the H2 molecule(s) carries(y) a partial charge, indicating that charge transfer takes place from the host Chemical">carbon atoms of the <Chemical">span class="Chemical">C50 cage to hydrogen atoms. Consequently, there exist electrostatic Coulomb forces between a H2 molecule(s) and C50.

The total Coulomb force between all H atoms of the H2 molecule(s) and all C atoms of the <Chemical">span class="Chemical">C50 fullerene cage was calculated using the following equationwhere f(nH2···C50) is the total force between the hydrogen molecule and all C atoms of the C50 cage and ε0 is the electric permittivity, and the sum is made over distinct pairs of particles (H, C) with charges (q, q) separated by a distance r.

Referring to the results shown in Table and based on the calculations using natural charges, the sum of the total attractive and repulsive forces between a H2 molecule(s) and the Chemical">C50 cage is nearly zero for H2@<Chemical">span class="Chemical">C50 (except for PBEPBE and to a lesser extent SVWN5) and non-negligible for 2H2@C50. Overall, these results indicate that only one H2 molecule encapsulated inside C50 is almost unaffected by the fullerene cage. However, for 2H2@C50 where 2H2 molecules are encapsulated inside C50, the attractive forces between H2 molecules and all C atoms of the fullerene cage lead to deformation of the fullerene cage. This is reflected by insignificant variation of the distance in the H2 molecule and particularly in the C50 cage in H2@C50, while for 2H2@C50, there clearly exists some variation of C–C bond lengths, where some are enlarged, whereas others are reduced, leading to an overall expansion (although slight) of the fullerene cage (see Figure and Table ).

On the other hand, the total Coulomb energies that are calculated based on the forces between H2 and Chemical">C50, although negligible, may give an idea on how stronger the interaction between H2 and the <Chemical">span class="Chemical">C50 cage is. Indeed, all calculations show that the total Coulomb energy is attractive in both complexes H2@C50 and 2H2@C50. Moreover, the Coulomb energy in 2H2@C50 is relatively much stronger than that in H2@C50 (see Table ), which further strengthens the fact that insertion of the 2H2 molecule will result in an expansion of the fullerene cage, while only one H2 does not affect.

Strain Energies

An exciting approach to understand the mechanism of host–guest interactions in H2@Chemical">C50 and <Chemical">span class="Chemical">2H2@C50 is to calculate their SEs. By definition, the energy stored in a body because of deformation is called the SE.[77] The SEs stored in nH2@C50 (n = 1 and 2) systems are determined by the overall summation of the SE of the host (C50) and that of the guest (nH2), as shown in eq

It is therefore expected from eq that the value of SE will be positive because both (E(Chemical">C50)Chemical">sp–complex – E(<Chemical">span class="Chemical">C50)opt) and (E(nH2)sp–complex – nE(H2)opt) terms are, in principle, positive.

On the other hand, eq can be written aswhere E(Chemical">C50)Chemical">sp–complex and E(nH2)Chemical">sp–complex correChemical">spond to the single-point energy calculations of the <Chemical">span class="Chemical">C50 cage and nH2 of the fully optimized endohedral complex, respectively, while E(C50)opt and E(H2)opt are the fully optimized energies of the isolated fullerene C50 and H2 molecules, respectively.

The addition of E(nH2@Chemical">C50)opt in both si<Chemical">span class="Chemical">des of eq leads towhere the first part is the complexation energy ΔE calculated as in eq (see the data listed in Table ), while the second part is referred to as ΔEsp.

Hence, the SE can be written as

Referring to eq , it is legitimate to expect that ΔEsp will be more negative than ΔE, in order to obtain a positive value for SE.

ΔE, ΔEsp, and SE were computed, and the results are listed in Table .

Table 7

SVWN5, PBEPBE, B3LYP, B2PLYP, B2PLYPD, and MP2 Complexation Energies (ΔE and ΔEsp) and SE of nH2@C50 Calculated Using the 6-311G(d,p) Basis Seta

	H2@C50			2H2@C50
	ΔE	ΔEsp	SE	ΔE	ΔEsp	SE
SVWN5	–12.26	–8.96	–3.30	3.43	–2.24	5.67
PBEPBE	–1.73	0.80	–2.53	26.83	18.48	8.35
B3LYP	4.46	4.22	0.24	35.37	25.96	9.41
B2PLYP	–5.47	–1.13	–4.34	24.04	15.07	8.97
B2PLYPD	–11.52	–7.24	–4.28	9.74	1.23	8.51
MP2(fc)	–8.03	–8.04	0.01	5.74	–2.67	8.41

All values are in kcal mol–1.

Surprisingly, not all calculated SEs are positive for H2@Chemical">C50 and not all of the ΔEChemical">sp values are more negative than their reChemical">spective ΔEs, as can be clearly noticed from Table . This “wrong” trend has also been found for several methods, when calculating ΔE, that is, finding positive values instead of negative ones, as can be seen in Table and by referring to eq . Because the calculations have been performed at 0 K and the host, guest, and complex all have Chemical">spin multiplicities equal to 1, the only explanation is that the method giving a positive value for ΔE is not adequate for this kind of system. Similarly, ΔEChemical">sp has to be, in principle, more negative, but what was found was the opposite. This means that the method(s) used is(are) again not adequate. Therefore, the negative values obtained in Table for SE and the wrong trend of ΔEChemical">sp with reChemical">spect to ΔE can be explained by nonadequacy of the method used. As can be seen in Table , the only exception was MP2, where ΔEChemical">sp is slightly more negative than ΔE and SE in slightly positive, proving that this method is appropriate for this kind of system and that the H2@<Chemical">span class="Chemical">C50 complex may be a metastable complex.

As for the Chemical">2H2@<Chemical">span class="Chemical">C50 complex, the overall trends in SE further confirm that it is thermodynamically unstable, in line with the obtained positive complexation energies shown in Table .

Dispersion Energies

When the interaction of nonpolar molecules is considered, dispersion is usually the dominant attractive force. Dispersion is a result of interactions between instantaneous charge fluctuations. The leading contribution to dispersion is due to an instantaneous dipole–dipole interaction, which gives rise to 1/r6 dependence of the DE. This type of intermolecular interaction could be considered in understanding the mechanism of host–guest interactions in H2@Chemical">C50 and <Chemical">span class="Chemical">2H2@C50. The DE of nH2@C50 (n = 1 and 2) systems is determined by the difference between the calculated complexation energies at empirical dispersion-corrected variation (B2PLYPD and B2PLYP), as shown in eq where DE, ΔEB2PLYPD, and ΔEB2PLYP correspond to the DE of the fully optimized nH2@C50 complex, calculated B2PLYPD complexation energy, and calculated B2PLYP complexation energy of nH2@C50 (n = 1 and 2), respectively. The DEs for H2@C50 and 2H2@C50 structures are given in Table .

Table 8

Computed DEsa for H2@C50 and 2H2@C50 Calculated Using the 6-311G(d,p) Basis Set with and without the BSSE Correction

	ΔEB2PLYP	ΔEB2PLYPD	DE
H2@C50	–3.05	–9.08	–6.03
2H2@C50	29.63	15.35	–14.28

All values are in kcal mol–1.

Table shows that the DE increases with increasing the number of Chemical">hydrogen molecules inside the <Chemical">span class="Chemical">fullerene cage. This is reflected in a significant decrease in the complexation energy of H2@C50 and 2H2@C50 systems, particularly the 2H2@C50 one. Although the value of the DE in the 2H2@C50 complex is rather higher (DE = −14.28 kcal mol–1), it was insufficient to bring the complex into a thermodynamically stable state.

On the other hand, it is well known that MP2 and MP4 levels of Møller–Plesset perturbation theory do account for DE.[78] It was also reported that MP2 overestimates the DE.[78] Therefore, the good agreement between MP2(fc)/6-311G(d,p) and B2BLYP/6-11G(d,p) values of −2.89 and −3.04 kcal mol–1, respectively, indicates that the true contribution of dispersion in H2@Chemical">C50 is relatively small, while the empirical diChemical">spersion correction in B2PLYPD somewhat overestimates the attractive diChemical">spersion.

Bond Dissociation Energy

For calculation consistency, the structure of Chemical">C50 (having <Chemical">span class="Chemical">D5h symmetry) that has been optimized at the MP2(fc)/6-311G(d) level was used in all calculations considered in this section. The expansion was made symmetrically, that is, following a breathing mode (all the origin-to-carbon distances are elongated or shortened with the same increment).

For computation limitations, the Chemical">BDE of <Chemical">span class="Chemical">C50 was calculated at PBEPBE and SVWN5 levels of DFT using the 6-311G(d) basis set, because both of them predicted similar trends as MP2 with regard to the stability of complexes H2@C50 and 2H2@C50, and also at B3LYP/6-311G(d) for comparative purpose. Here, the BDE is defined as the enthalpy change in the gas phase of C50 into 50 carbon atoms (C50 → 50 × C) calculated at a temperature of 298.15 K and a pressure of 1 atm.

The eChemical">nthalpy of formation for <Chemical">span class="Chemical">C50 fullerene was calculated using the following equationwhere ZPE is the zero-point energy correction; Htrans., Hrot., and Hvib. are the standard temperature correction terms calculated using the equilibrium statistical mechanics with a harmonic oscillator and rigid rotor approximations; RT is the ability to do pressure–volume work (PV).[76]

The total Chemical">BDE of <Chemical">span class="Chemical">C50 is defined as

Moreover, the Chemical">BDE per bond is given bywhere <Chemical">span class="Chemical">BDE is the total bond dissociation energy, BDE/bond is the BDE per bond, HC is the enthalpy of formation of the optimized C50, and HC is the enthalpy of formation of the carbon atoms. The number of C–C bonds is represented by k (which is equal to 75 bonds in C50).

Figure shows the average Chemical">BDE per bond (<Chemical">span class="Chemical">BDE/bond) calculated at PBEPBE/6-311G(d,p), SVWN5/6-311G(d,p), and B3LYP/6-311G(d) levels of DFT as a function of expansion (ε).

Figure 3

Average C–C BDE per bond (BDE/bond) of C50 calculated at PBEPBE/6-311G(d), SVWN5/6-311G(d), and B3LYP/6-311G(d) levels of DFT as a function of expansion (ε).

Average C–C Chemical">BDE per bond (<Chemical">span class="Chemical">BDE/bond) of C50 calculated at PBEPBE/6-311G(d), SVWN5/6-311G(d), and B3LYP/6-311G(d) levels of DFT as a function of expansion (ε).

According to Figure , the calculated Chemical">carbon–<Chemical">span class="Chemical">carbon BDE of C50 is 111.5 kcal mol–1 at the PBEPBE level (127.1 kcal mol–1 at SVWN5 and 102.6 kcal mol–1 at B3LYP), which is larger than that of ethane =  and smaller than that of ethylene .[79] This finding is consistent with the average C–C bond length of C50 (see r(average) in Tables –4) that is shorter than that of carbon–carbon single-bond distance (C–C = 1.54 Å) and longer than that of carbon–carbon double-bond distance (C=C = 1.34 Å).[80] According to Tables –4, the average C–C bond length is of the order of 1.43–1.44 Å, regardless of whether C50 is isolated or the H2 molecule(s) is(are) encapsulated inside C50. Furthermore, none of the C–C bond length (Rmin or Rmax) reaches the limit of single or double bonds. The high value of the BDE for C50 accounts for the high stability of the C50 molecule. Consequently, considering the high stability of C50, the energy required for the bond elongation should be prohibitive. For example, the calculated energy required for C50 expansion by only ε = 6.2% is 390.4 kcal mol–1 at the PBEPBE/6-311G(d) level (359.8 kcal mol–1 at B3LYP and 298.9 kcal mol–1 at SVWN5). On the other hand, ΔH = Hexpand – Hrelaxed, where ΔH, Hexpand, and Hrelaxed are the energy required for C50 expansion and the enthalpy of formation of the expanded and relaxed C50 fullerene, respectively. Therefore, the mechanism proposed for C50 expansion could be realized only in a hypothetical system. In a real situation, as shown in Tables –4, the nH2@C50 (n = 1 and 2) complex is characterized by the same carbon–carbon bond lengths of a free C50 cage optimized within the same method of calculation.

Conclusions

In this work, we reported DFT (pure and hybrid), HF, MP2, MP4, B2PLYP, and B2PLYPD investigation on the possibility of encapsulation of H2 and Chemical">2H2 molecules inside the <Chemical">span class="Chemical">C50 fullerene cage. Based on MP2-, MP4-, B2PLYP-, and B2PLYPD-calculated complexation energies, the stability of the H2@C50 complex is predicted, which was featured by negative complexation values. SVWN5 also predicted a stable H2@C50 complex but with overestimation, while PBEPBE gave slight stability. Taking into consideration the BSSE correction, it was found that the DFT functionals and long-range-corrected functionals often tend to underestimate the stabilization energy in the H2@C50 complex, while the MP2 method overestimates it. However, for the complex with two encapsulated hydrogen molecules, 2H2@C50, all methods used gave positive complexation energies, proving the instability of this complex. Among the various methods used in this research, MP2 and B2PLYPD, in conjunction with the triple-zeta Pople-style 6-311G(d,p) basis set, provide the most reliable results in predicting the stability of nH2@C50 complexes.

These calculations were further supported by consideration of host–guest interactions such as charge transfer, Coulomb forces, Coulomb energy, SE, DE, and Chemical">BDE of the <Chemical">span class="Chemical">C50 cage. Indeed, by assuming that one H2 molecule is to be encapsulated with the equilibrium H–H internuclear distance and by summing up this distance with the van der Waals distances of C and H atoms, it was found that the C50 fullerene cage can merely accommodate one H2 molecule. However, even if we consider that H2···H2 distance is that of the optimized complex, the encapsulation of the second H2 molecule is thermodynamically an unfavorable process. These findings are proven by MP2, MP4, B2PLYP, and B2PLYPD calculations.

The analysis of the host–guest interaction energies reveals that the main stabilizing effects are always due to the SE and DE, while the electrostatic Coulomb energy contributes negligibly to the stabilization of the complexes.

On the other hand, the calculated Chemical">carbon–<Chemical">span class="Chemical">carbon BDE of C50, which is larger than that of ethane and smaller than that of ethylene, is consistent with the average C–C bond length of C50 that is shorter than that of carbon–carbon single-bond distance and longer than that of carbon–carbon double-bond distance. The high value of the BDE for C50 accounts for the high stability of the C50 molecule. Consequently, the nH2@C50 (n = 1 and 2) complex is characterized by the same carbon–carbon bond lengths of an isolated C50 cage.

Finally, it is apparent from this work that we also did a benchmark study for the encapsulation of the H2 molecule inside the Chemical">C50 fullerene cage by employing several methods and basis sets. With no doubt, the larger basis sets significantly improve the quality of the results and can even change the trends (case of PBEPBE). However, only methods which involve treatments of diChemical">spersion interactions between the host and guest are found to be adequate for encapsulation of H2 inside <Chemical">span class="Chemical">fullerene cages.