by Hydrogen Clusters at

Subscriber access provided by UNIV OF DURHAM

Article 2

60

Physisorption of H on Fullerenes and the Solvation of C by Hydrogen Clusters at Finite Temperature: A Theoretical Assessment Florent Calvo, Ersin Yurtsever, and Adem Tekin J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.8b00163 • Publication Date (Web): 16 Feb 2018 Downloaded from http://pubs.acs.org on February 22, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Physisorption of H2 on Fullerenes and the Solvation of C60 by Hydrogen Clusters at Finite Temperature: A Theoretical Assessment F. Calvo,∗,† E. Yurtsever,‡ and A. Tekin¶ LiPhy, Université Grenoble 1 and CNRS UMR 5588, 140 Avenue de la Physique, 38402 St Martin d’Hères, France, Koç University, Chemistry Department, Rumeli Feneri Yolu, Sariyer 34450, Istanbul, Turkey, and Informatics Institute, Istanbul Technical University, 34469 Maslak, Istanbul, Turkey E-mail: [email protected]

∗

To whom correspondence should be addressed LiPhy, Université Grenoble 1 and CNRS UMR 5588, 140 Avenue de la Physique, 38402 St Martin d’Hères, France ‡ Koç University, Chemistry Department, Rumeli Feneri Yolu, Sariyer 34450, Istanbul, Turkey ¶ Informatics Institute, Istanbul Technical University, 34469 Maslak, Istanbul, Turkey †

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract The interaction between hydrogen and carbonaceous nanostructures is of fundamental interest in various areas of physical chemistry. In this contribution we have revisited the physisorption of hydrogen molecules and H2 clusters on fullerenes, following a first-principles approach in which the interaction is quantitatively evaluated for the C20 system using highlevel electronic structure methods. Relative to coupled cluster data at the level of single, double, and perturbative triple excitations [CCSD(T)] taken as benchmark, the results for rotationally averaged physisorbed H2 show a good performance of MP2 variants and symmetryadapted perturbation theory (SAPT), but significant deviations and basis set convergence issues are found for dispersion-corrected density functional theory (DFT). These electronic structure data are fitted to produce effective coarse-grained potentials for use in larger systems such as C60 –H2 . Using path-integral molecular dynamics, the potentials are also applied to parahydrogen clusters solvated around fullerenes, accross the regime where the first solvation shell becomes complete and as a function of increasing temperature. For C60 our findings indicate a sensible dependence of the critical solvation size on the underlying potential. As temperature is increased, a competition is found between the surface and radial expansions of the solvation shell, one molecule popping away at intermediate temperatures but getting reinserted at even higher temperatures.

1 Introduction The interaction between molecular hydrogen and carbon surfaces has attracted a lot of attention in the context of hydrogen storage, with the hope that nanoporous carbon structures could hold significant quantities of H2 and meet the storage density requirements for industrial applications. 1–5 From a more fundamental perspective, hydrogen in contact with carbonaceous grains has been the subject of intense research owing to the possibility that H2 might be formed by surface recombination of atomic hydrogen. 6–8 Hydrogen weakly binds to clean (unreactive) carbon surfaces, and the accurate determination

2


Page 2 of 28

Page 3 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of the binding energy from first principles becomes difficult if the surface is extended or disordered. Most treatments in molecular simulation have used empirical force fields adjusted on experimental properties, with the aim of conducting possibly large scale investigations on extended media such as nanotubes. 9,10 Hydrogen adsorption on graphite has been studied theoretically using ab initio approaches 15 sometimes in combination with other methods such as kinetic Monte Carlo, 16 again in relevance with hydrogen energy or astrophysical motivations. Exohedral adsorption of H2 on pristine and doped fullerenes has been considered by Wang and coworkers, 17,18 however these authors did not account for dispersion forces in their density functional theory (DFT) calculations. Endohedral hydrogen fullerene complexes have also been addressed along similar lines, as well as semi-empirical 11 or quantum chemical 12,13 calculations, notably to investigate the transition from physisorption to chemisorption as hydrogen density increases inside the fullerene. 14 In the limit of individual molecules, highly accurate potential energy surfaces have been designed to address the quantum dynamics of H2 impinging on graphite surfaces 19 or inside fullerene cages. 20 Recently, mass spectrometry experiments under the extreme cryogenic conditions provided by helium nanodroplets were performed with a controlled number of hydrogen molecules coating the + 21 fullerene cations C+ The variations of mass abundances revealed shell completions at 60 and C70 .

49 and 51 molecules for these two molecules and for both H2 and D2 coatings. These measurements could be partly rationalized using DFT calculations performed for individual molecules, 21 and more systematic path-integral molecular dynamics (PIMD) simulations for hydrogen molecular clusters based on a coarse-grained potential for H2 derived from electronic structure calculations. 22 Density-functional calculations by Kaiser and coworkers 21 reported some significant variations of the H2 -C60 interaction depending on the functional, reaching about 50% in magnitude for the binding energy between PBE0 and wB97xD. However, in larger clusters the geometrical packing schemes constrained by the H2 –H2 interaction were found to produce a relatively minor influence of such variations in the H2 -C60 physisorption energy on the number of hydrogen molecules that the fullerene can accomodate in its first solvation shell. 22 The apparently poor ability of DFT to determine accurate binding energies for such non-

3



covalent complexes is not entirely surprising, even considering that dispersion forces were in principle accounted for in the two aforementioned functionals. The goal of the present contribution is to discuss alternative electronic structure methods that could be more appropriate for the problem of hydrogen-fullerene exohedral interaction, and use them to predict and compare the solvation shell size in (neutral) C60 . Because buckminsterfullerene has too many electrons, we have considered the smaller dodecahedral cage of C20 as our reference system, enabling coupled cluster calculations at the level of single, double, and perturbative triple excitations. Besides CCSD(T), our investigation includes second order perturbation theory (MP2) and its spin-component scaled variants, symmetry-adapted perturbation theory (SAPT) based on DFT, as well as DFT itself that was revisited to address the specific behavior toward basis set convergence, here employing the wB97xD functional as in Ref. 21. The electronic structure calculations cannot be straightforwardly scaled to larger fullerenes or assemblies of hydrogen molecules, and we have transformed these data into more practical ab initio potentials, keeping a coarse-grained (structureless) character for the hydrogen molecules that is suitable for gas phase systems. The various potentials generated from the different electronic structure methods were employed to simulate hydrogen clusters around C20 and C60 at finite temperature, extending our earlier work 22 to the determination of possible phase changes in hydrogen clusters in presence of such dopants. Unsurprisingly, variations in the interaction energies inferred from electronic structure methods are reflected on the collective properties of hydrogen clusters, although no qualitatively different behavior was detected. The article is organized as follows. In the next section we detail the different electronic structure methods used in this work for C20 -H2 , the semiempirical potential employed to treat larger systems, and we provide some brief information about the PIMD methodology. The results are presented and discussed in Sec. 3, before some concluding remarks are finally given in Sec. 4.

4


Page 4 of 28

Page 5 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2 Methods 2.1

Electronic structure calculations

The intermolecular interaction energy can be calculated with a supermolecular approach at any level of theory. CCSD(T) is expected to be the most accurate applicable method for the calculation of interaction energies for the present compounds. However, its extreme computational demand forces us to find alternative methods which are much cheaper but hopefully nearly as accurate as CCSD(T). The second order Møller-Plesset perturbation theory (MP2), spin-component scaled MP2 (SCS-MP2) 31 and SCS-MI-MP2 33 or dispersion-corrected density functional theory (DFT-D) are currently among the best alternatives to CCSD(T). In supermolecular interaction energy calculations, the computed energy is just a single number and therefore it does not contain any information about the nature of the interaction. Symmetryadapted perturbation theory (SAPT), 34,35 an alternative method to the supermolecular approach, can also be used for the calculation of interaction energies. In SAPT, the total interaction en(1)

ergy is decomposed into physically distinct terms such as first-order electrostatic (Eel ), second(2)

(2)

(1)

(2)

order induction (Eind ) and dispersion (Edisp ) and their exchange counterparts (Eexch , Eexch−ind and (2)

Eexch−disp ). Effects of higher order contributions are incorporated with a δ(HF) term, which is the difference between the supermolecular Hartree-Fock (HF) interaction energy and the electrostatic, induction and their exchange counterparts energies obtained from HF monomer density matrices and HF response functions. Recently, a SAPT implementation using DFT description of the monomer properties and refered to as DFT-SAPT was developed by Hasselmann and Jansen, 36–39 and was found to show superior scaling bahvior compared to the original many-body version of SAPT. 34 This implementation was further equipped with the density-fitting approximation called as DF-DFT-SAPT. 40 The MOLPRO quantum chemistry program package 41 was employed for the supermolecular (HF, MP2, SCS-MP2, SCS-MI-MP2 and CCSD(T)) and DF-SDFT-SAPT interaction energy calculations on the C20 –H2 potential energy surface. Dunning’s aug-cc-pVXZ (aVXZ) basis sets with

5



Page 6 of 28

X=D and T were employed in the interaction energy computations. Weigends’s cc-pV(X+1)Z JK-fitting basis 42 in conjunction with an aug-cc-pVXZ basis set was used in DF-HF and DF-DFTSAPT. The corresponding aug-cc-pVXZ MP2-fitting basis set 43 was employed in DF-MP2, DFSCS-MP2 and DF-SCS-MI-MP2. Asymptotically corrected PBE0AC 36 and LPBE0AC 40 exchangecorrelation (xc) potentials in combination with the adiabatic local density approximation (ALDA) xc-kernel were utilized in DF-DFT-SAPT computations. All supermolecular interaction energies were corrected for basis set superposition errors (BSSE) with the help of the counterpoise correction (CP). 44 DF-DFT-SAPT interaction energies are free from the BSSE. Intermolecular interaction energies computed with aVDZ and aVTZ basis sets were extrapolated to the complete basis set (cbs) limit using a two-point formula proposed by Bak et al. 45 and described elsewhere. 46 These extrapolated energies were then employed in the fitting process for the generation of the empirical potential energy functions.

2.2

Semiempirical modeling

The explicit potential developed to model H2 molecules in contact with polycyclic aromatic hydrocarbons 23 and fullerenes 22 was naturally employed here to describe assemblies of hydrogen molecules around C20 and C60 . It relies on the structureless approximation for H2 , which is treated as a pointlike particle owing to the free rotation it experiences in its J = 0 state at the low pressures relevant for isolated systems. Different potentials are available to model the interaction between hydrogen molecules, 24,25 and the widespread Silvera-Goldman (SG) potential 25 was chosen here in its ’isolated’ parametrization. The interaction V (R) between each H2 molecule and the carbonaceous dopant is treated assuming a pairwise additive form over all carbon atoms, with a similar functional form as in the SG model. Denoting here R as the configuration of the system with atomistic details for the fullerene, the potential felt by hydrogen molecule i thus reads

Vi =

X

[Vrep (rij ) + Vdisp (rij )]

j∈CN

6


(1)

Page 7 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


with the repulsive and dispersion terms written as

Vrep (rij ) = A exp(−brij ) "

Vatt (r) = −fcut (rij ) 

(2)

C6 C8 + 8 6 rij rij

#

(3)

!2   , rij ≤ rc ,

rc fcut (rij ) = exp − −1 rij = 1,

rij > rc

(4)

This simple potential only contains four parameters (A, b, C6 and C8 ), in addition to the cutoff parameter rc . Since we are only dealing with neutral fullerenes, the polarization contribution considered in our earlier work for cationic dopants 22,23 was neglected here. In particular, we did not consider more sophisticated descriptions of the multipolar distribution in buckminsterfullerene, 26 although we anticipate that this would affect the results presented below only very marginally.

2.3 Path-integral molecular dynamics We do not give an extensive presentation of the PIMD method, which has been reviewed recently 27 notably in the context of gas phase systems, 28 but only provide its aspects relevant to the present systems. In short, PIMD is a rigorous computational method to incorporate the quantum mechanical nature of the nuclei, except for exchange statistics that we neglect here. In practice, the temperature is chosen sufficiently high to ensure that the bosonic nature of parahydrogen can be safely ignored, that is above 1 K. 29 The PIMD method replaces the classical particles by a set of P monomers or beads that are consecutively connected to each other in a cyclic fashion through appropriate harmonic bonds, forming a ring-polymer whose extension should be understood as describing the extent of vibrational delocalization. The method becomes exact in the limit where the Trotter number P becomes infinite, but for the present application a reasonable value of P = 128 was employed. In practice, PIMD trajectories proceed similarly as classical (thermostatted) molecular dynam-

7



Page 8 of 28

ics, a Nosé-Hoover thermostat being here attached to each of the P replicas of the system. The velocity Verlet algorithm is employed to integrate the equations of motion with a time step of 0.2 fs.

3 Results and discussion 3.1

C20 -H2 interaction from electronic structure calculations

The parameters of the C20 –H2 potential described in Eqn. (1-4) were fitted to the quantum chemical calculations employing different methodologies and basis sets. The single radial parameter of the interaction potential is chosen as the distance between the center of masses of H2 and C20 and varied from 4 to 9 Å. At selected points, H2 was placed along three orthogonal axes and the interaction energies were calculated separately while keeping C20 and H2 frozen. After correcting for possible basis set superposition errors, these energies were then averaged over the three directions to calculate the final potential energy. At these low temperatures, energy differences due to the rotational motion of H2 are negligible. The parameters of the semiempirical potentials were eventually fitted by an adaptive Monte Carlo method. Such obtained parameters as well as our original parameters 22 are given in Table 1. Table 1: Optimal parameters of the H2 -carbon potential of Eqn. (2–3) representing best the effective C20 –H2 interaction, adjusted on dedicated electronic structure calculations or borrowed from Ref. 22. All parameters are given in atomic units. Fitted data CCSD(T) DFT-SAPT MP2 wB97xD Ref. 22

A 9.492 16.853 18.540 48.540 73.828

b 1.577 1.685 1.686 1.919 1.999

C6 24.129 26.651 32.098 20.091 19.406

C8 382.534 221.969 277.001 96.186 139.526

Quantum mechanically calculated points and fitted potentials are represented in Figs. 1–3 for the C20 –H2 complex. Since the fittings were performed for one-dimensional functions, the Monte Carlo fitting performs very well for all cases. The most accurate potential function from CCSD(T) 8


Page 9 of 28

0

Energy (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


CCSD(T) DFT-SAPT MP2 wB97xD potential

-2

-4 4

5

6

7

8

9

Distance from center (Å) Figure 1: (Color Online) Potential energy curves for C20 –H2 after averaging on H2 orientations, using different electronic structure methods and the reference empirical potential. The continuous lines are the optimized potentials for each method using the same empirical form as given in the text. The results are shown for extrapolation to complete basis set limit, using using basis sets that generally rely on Dunning’s aug-cc-PVXZ basis sets, as described in the text. is used as the standard for testing the quality of the various other methods. In Fig. 1, we present the results from different methodologies: DFT with wB97xD which is one increasingly popular functional for long range interactions, MP2, DFT-SAPT, CCSD(T) and our reference empirical potential that we denote below simply as the reference potential. The minimum energy from CCSD(T) is −2.44 kJ/mol at an equilibrium distance of 4.9 Å. As expected, MP2 overestimates bonding in non-covalent systems 30 with a binding energy of −3.51 kJ/mol at 4.8 Å. In contrast, wB97xD predicts a slightly lower binding energy of −2.17 kJ/mol but approximately the same equilibrium distance as CCSD(T). DFT-SAPT as well as our potential agree with CCSD(T) results very well. For the long-range distances, only DFT-SAPT stays below CCSD(T) results.

9



1

0

Energy (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 28

-1

CCSD(T) MP2 SCS-MP2 SCS-MI-MP2

-2

-3

-4 4

5

6

7

8

9

Distance from center (Å) Figure 2: (Color Online) Potential energy curves for C20 –H2 after averaging on H2 orientations, using MP2 and related methods. Here the lines joining the symbols are guides to the eye. These calculations employed the aug-cc-PVXZ MP2-fitting basis set and were extrapolated to complete basis set limit (see text for details). In Fig. 2 we compare the improvements made to the MP2 methods. SCS-MP2 averages over the singlet and triplet energies with appropriate weight of 1.2 and 1/3. This method was shown to give much better interaction energies for non-covalent systems. However, it is known that this correction may be overcompensating the binding and it is recommended that the weights could be optimized for a given system. 31,32 Here, SCS-MP2 corrects the binding energy from −3.52 kJ/mol to −2.18 kJ/mol which is now lower than the CCSD(T) result. Upon rearranging the weights of the singlet and triplet states as described in section 2.1, a much better fit is now obtained for SCSMI-MP2 both in terms of the strength of the bonding and the equilibrium distance which become −2.35 kJ/mol and 4.9 Å, respectively. The effect of the basis set size is given in Fig. 3 specifically for DFT calculations with the 10


Page 11 of 28

1

0

Energy (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


-1

aug-cc-pVDZ aug-cc-pVTZ extrapolation

-2

-3

-4 4

5

6

7

8

9

Distance from center (Å) Figure 3: (Color Online) Potential energy curves for C20 –H2 after averaging on H2 orientations, using density functional theory with the wB97xD functional and different basis sets. wB97xD functional. As aforementioned, the basis sets aug-cc-pVDZ and aug-cc-pVTZ were used for a complete basis set extrapolation, keeping in mind that the dispersion correction of this functional is analytical thus basis set independent. Normally it is expected that, upon increasing the basis set, the binding energies would also increase. This is the behavior observed for all methods employed in the present work with the exception of this DFT/wB97xD approach. In this case, the extrapolated binding energy is about 1.5 kJ/mol lower than obtained with the small basis set augcc-pVDZ. Unfortunately, CCSD(T) calculations could not be performed with these larger basis sets so we cautiously do not conclude whether this is a peculiarity of the present DFT method or the C20 –H2 system. Dispersion-corrected DFT with the wB97xD method thus appears to entail some possible issues associated with basis set convergence. We performed additional calculations with the B3LYP-D3 11



method near the minimum geometry and did not find this behavior, but rather weak basis set effects. In addition to not being general to dispersion-corrected functionals, it is also possible that this unexpected behavior could be an artifact of the BSSE correction. This could then indicate that the couterpoise method fails for this specific system and method, but in any case this would require further investigation. 1

(a) PBE0AC

1

-1

(1)

Eel

Eexch

(1)

(2)

Eind

Eind-exch

-2

Edisp

(2)

5

6

-1

-2

(2)

Edisp-exch HF Eint

-3 4

(b) LPBE0AC

0

energy (kJ/mol)

0

energy (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 28

(2)

-3 4

7

distance (Å)

5

6

7

distance (Å)

Figure 4: (Color Online) DFT-SAPT energy contributions for the C20 –H2 orientation with H2 aligned with the intermolecular axis, using (a) PBE0AC and (b) LPBE0AC as the density func(2) (2) (1) (1) tionals. Eel and Eexch are shown with long dashed lines, Eind and Eind−exch with short dashed (2) (2) lines, Edisp and Eexch−disp with dotted lines, δ(HF) with dot-dashed lines, and Eint with thick solid lines. Figs. 4(a) and 4(b) show the energy components obtained with DFT-SAPT(PBE0AC) and DFTSAPT(LPBE0AC), respectively, using the aug-cc-pVTZ basis set as a function of intermolecular (1)

(2)

(2)

distance, assuming here the hydrogen molecule to be aligned along the axis. Eel , Eind and Edisp were always obtained as negative. δ(HF) was also negative except for the short intermolecular (1)

(2)

separations. The major repulsive contribution is coming from Eexch . Eind is significant especially (2)

at short R while Edisp is active at longer R. DFT-SAPT(LPBE0AC) lowers the total interaction energy compared to DFT-SAPT(PBE0AC) by almost 0.5 kJ/mol. This is mainly due to the more (2)

pronounced Edisp term in DFT-SAPT(LPBE0AC).

12


Page 13 of 28

3.2

Application to C60 -H2 1

0

Energy (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


-1

-2

-3

CCSD(T) DFT-SAPT MP2 wB97xD potential

-4

-5 5

6

7

8

9

10

Distance from center (Å) Figure 5: (Color Online) Potential energy curves for C60 –H2 predicted by the different potentials fitted to reproduce the various reference electronic structure data for C20 –H2 . Using the semiempirical potential energy functions fitted to reproduce the quantum chemistry results, we calculated the interaction between C60 and H2 . Buckminsterfullerene C60 is composed of 32 faces where 20 are hexagons and 12 are pentagons, hexagonal sites bind slightly more than pentagonal sites. Fig. 5 shows the binding energy for H2 approaching C60 along a sixfold axis, again assuming rotational averaging for the H2 axis. The same trends as in C20 –H2 case are found. The MP2 curve has a much stronger binding energy and a slightly shorter equilibrium distance, while wB97xD yields the smallest binding energy. The other three methods display similar characteristics. The equilibrium distance shifts from 4.7 Å for C20 to to 6.3 Å for C60 , consistently with the larger radius of this fullerene (1.55 Å higher than that of C20 ). In a similar fashion, the potential goes to zero at a much slower rate and even at 10 Å, the binding energy is still around 13



Page 14 of 28

−0.2 kJ/mol. Using the above potentials, we have also considered the different adsorption sites for an hydrogen molecule on C60 , namely hexagonal and pentagonal facets, as well as the effect of zero-point energies (ZPE) for such light adsorbant species. Table 2 shows the physisorption energies for the five parametrizations of the potential and the two sites, with and without ZPE corrections in the harmonic approximation. Table 2: Physisorption energies of H2 on hexagonal and pentagonal sites of C60 , as predicted by the empirical potential fitted to reproduce different electronic structure calculations, in the classical case and with harmonic zero-point energy (ZPE) corrections. All energies are in kJ/mol. Fitted data CCSD(T) DFT-SAPT MP2 wB97xD Ref. 22

Hexagonal site -3.154 -3.592 -4.887 -2.931 -3.260

Pentagonal Hexagonal site site + ZPE -2.869 -2.080 -3.257 -2.373 -4.421 -3.443 -2.642 -1.713 -2.927 -1.915

Pentagonal site + ZPE -1.913 -2.169 -3.134 -1.554 -1.726

All potentials consistently predict H2 to bind more strongly to the hexagonal site by about 10%, and the same ordering between underlying ab initio methods is of course found as in Fig. 5. The inclusion of zero-point correction does not change this ordering, nor the magnitude of the greater energetic stability of the hexagonal site. However, vibrational delocalization is found to be quantitatively significant, reducing the physisorption energy by more than 30% for all methods. This confirms the importance of describing nuclear motion on a quantum mechanical footing when treating assemblies of hydrogen molecules around fullerenes.

3.3

Stepwise solvation of fullerenes by parahydrogen

The solvation of C20 and C60 has been studied by PIMD calculations at 2 K with stepwise increments in the number of hydrogen molecules coating the cluster, across the range where the first solvation shell becomes complete. For C20 (H2 )n , n was chosen between 20 and 40 while the range 14


Page 15 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


for C60 (H2 )n was taken as 40–60, in relation with mass spectrometry experiments that suggested shell completion near n = 50. 21 Based on our previous experiences with PIMD of hydrogen clusters with carbonaceous dopants 22,23 we choose the Trotter number P = 128 and a time step of 0.2 fs; each structure was equilibrated by 2 × 105 time steps and the averages accumulated over the subsequent 106 time steps. Several structural and dynamical properties were considered for analysis, namely the shell size, shell radius, and Lindemann index. With these computational parameters we have tested all 5 sets of potential parameters reported in Table 1. At the low temperature considered, H2 clusters form a well-defined first solvation shell for both fullerenes. In other hydrocarbon molecules such as polycyclic aromatic hydrocarbons, 22 the second solvation shell is also distinguishable in the radial distribution functions and even the third shells starts to form provided that a sufficient number (typically n > 100) of H2 molecules are available. In case of the fullerenes, and owing to their near spherical symmetry, up to 4 distinct shells could be detected with 500 H2 molecules. 23 In the present study, we restrict ourselves to only the properties of the first shell and in the above mentioned ranges, where the shell is complete. The two fullerene dopants having all carbon atoms on a same sphere, two structural parameters can be calculated from the radial distribution functions. Here we considered the average distance and, by integrating only the first peak, the number of particles within the first shell. In Fig. 6(a,b), we represented both quantities for C20 (H2 )n as functions of the size n. For this dopant, hydrogen molecules initially bind prefentially near the centers of the 12 pentagonal faces, but the surface density remains low enough for the system to accomodate more molecules in its first shell, which is clearly seen in Fig. 6(a). The addition of molecules initially proceeds with a slight expansion in the shell radius [see Fig. 6(b)], until the shell is complete at n = 32. At this size, the number of molecules in the shell becomes stationary and the slope in the radial expansions exhibits a strong increase marking the appearance of an outer shell. These results are qualitatively independent of the potential used for the PIMD calculations, except for some differences in the shell radius that are consistent with the equilibrium distance between C20 and H2 in Fig. 1. MP2, in particular, systematically predicts significantly smaller solvation shells, while wB97xD slightly overshoots.

15



35

Shell size

(a) 30 CCSD(T) DFT-SAPT MP2 wB97xD potential

25

20

Average radius (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 28

6

(b)

5.8 5.6 5.4 5.2 20

25

30

35

Number of hydrogen molecules

40

Figure 6: (Color Online) Structural properties of C20 (H2 )n clusters versus size n, in the range where the solvation shell becomes complete, as predicted by PIMD simulations at 2 K using different parametrizations of the potential. (a) Shell size; (b) average radius. As seen in Fig. 7(a,b), in the case of C60 (H2 )n the changes in the average radii are essentially similar to those in C20 (H2 )n and as expected the smallest and largest radii are obtained from calculations with the MP2 and wB97xD potentials, respectively. Again there are two different slopes but the transition range is not as well defined as in C20 . At small sizes, hydrogen binds preferentially to the hexagonal and pentagonal sites, but once all of them (32) are filled there is still room for more molecules to fill the shell. For this dopant the shell size monotonically increases up to n = 51, but the significant fluctuations yield only a transition range 51–54 depending on size and the potential [see Fig. 7(b)]. Consistently with the relative strengths of the various potentials the shell has the shortest radius but the largest size for the MP2 parameters, opposite trends being found for the wB97xD parameters. The fluctuations in the shell size are partly due to the potential itself, but are also intrinsic size effects that were noted earlier. 23 In addition, they are also strengthened by temperature fluctuations, which we have explored by repeating the PIMD simulations in the range 2–15 K. For these simulations, the Trotter number 16


Page 17 of 28

55

Shell size

(a) 50 CCSD(T) DFT-SAPT MP2 wB97xD potential

45

40 7.2

Average radius (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(b) 7 6.8 6.6 6.4 40

45

50

55

Number of hydrogen molecules

60

Figure 7: (Color Online) Structural properties of C60 (H2 )n clusters versus size n, in the range where the solvation shell becomes complete, as predicted by PIMD simulations at 2 K using different parametrizations of the potential. (a) Shell size; (b) average radius. P was varied by keeping P T = 256 K as a constant. Besides the shell size and its radius, the Lindemann index was also considered as a measure of fluxionality in the system dynamics. This quantity was obtained from the centroid positions ¯ri of the hydrogen particles as h

2 rij i − h¯ rij i2 X h¯ 2 δ= n(n − 1) i 0.1). From a quantitative perspective, a similar phenomenology is obtained for the MP2 set of parameters, but the transition is delayed to 10 K and is accompanied with a much shorter radius. For comparison, we also show in Fig. 8 the results of classical simulations (corresponding to P = 1), performed under the same conditions for the specific set of parameters of Ref. 22 The differences with path-integral MD are striking, with the Lindemann index remaining below 10% throughout the entire temperature range T < 15 K, the shell size remaining essentially constant and equal to the system size, and a much lower average radius consistent here with the absence of radial expansion caused by vibrational delocalization. Melting below 10 K for this system thus 18


Page 19 of 28

partly originates from quantum mechanical effects. The trends in C60 (H2 )50 are generally parallel to those of C20 (H2 )32 (Fig. 9), except for one interesting difference. Again the MP2 potential stands out as producing similar variations in the Lindemann index [Fig. 9(a)], shell size [Fig. 9(b)] and average radius [Fig. 9(c)] as the other potentials, only delayed to higher temperatures. The striking features for this cluster is the nonmonotonic variations in the shell size and average radius above 7–8 K, with a dip in the shell size at 5–8 K depending on the potential, the Lindemann index showing significant increases only at the upper

Lindemann index

limit of this range. Equilibrium densities for this system are represented in Fig. 10(a-d) at various 0.25 (a) 0.2 0.15 0.1 0.05 0

Shell size

50 49.8

(b)

CCSD(T) DFT-SAPT MP2 wB97xD potential (classical)

49.6 49.4 49.2 49

Average radius (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


7 6.9

(c)

6.8 6.7 6.6 6.5 0

5

10

Temperature (K)

15

Figure 9: (Color Online) Structural properties of the C60 (H2 )50 cluster as a function of temperature, as predicted by PIMD simulations at using different parametrizations of the potential. (a) Lindemann index between centroids; (b) shell size; (c) average radius. Classical MD results using the potential parameters of Ref 22 are superimposed with orange stars . temperatures, with the classical structure shown in Fig. 10(a) and the mostly localized classical trajectory at 2 K illustrated in Fig. 10(b). The densities at 6 K [Fig. 10(c)] and 15 K [Fig. 10(d)] reveal that the dip is a manifestation of one molecule popping out of the shell (and only one) but 19



inserting itself back at higher temperatures once the underlying shell has expanded enough. This nonmonotonic size effect is robust against all parameter sets apart from some quantitative differences in its onsets, and is a manifestation of stronger finite size effects for this larger dopant, which is also supported by more larger fluctuations in the experimental mass abundances for the cationic 21 impurity C+ 60 .

(a) Classical GM

(b) Classical (2 K)

(c) PIMD (6 K)

(d) PIMD (15 K)

Figure 10: (Color Online) Classical and quantum structures of C60 (H2 )50 at various temperatures as predicted with the reference potential. The more complete isomerization of the hydrogen shell occurs once the floating molecule has reinserted into it, in the 8–10 K range depending on the potential. This higher melting temperature compared to C20 (H2 )32 is consistent with the stronger effective binding of hydrogen molecules to the larger 60-atom dopant. In the classical limit, additional simulations do not find evidence for such mechanisms, and the hydrogen shell remains solid like and locked to the fullerene throughout the entire temperature, without any molecule leaving the main shell. This behavior is very similar to helium clusters adsorbed on C60 , for which except near 32 atoms coverage the shell behaves solidlike when nuclear 20


Page 20 of 28

Page 21 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


motion is treated classically while a quantum mechanical treatment predicts a liquid. 48

4 Conclusions Exohedral absorption of hydrogen on carbon nanostructures is a very weakly energetic process in which the van der Waals complexes are only stable under cryogenic conditions below a few tens of Kelvin in the gas phase. The present computational work was aimed at revisiting this interaction between H2 molecules and small fullerenes, paying a particular attention at C20 for which a variety of quantum chemical methods can be applied. Using such a broad arsenal, we found SAPT-DFT and spin-scaled MP2 methods to perform best against coupled cluster benchmark calculations. In contrast, standard (unscaled) MP2 was found to overestimate the binding strength significantly while dispersion-corrected DFT appears to lack reliability with basis set convergence issues that prevent a clear extrapolation to complete basis set limit. An effective pairwise potential between structureless H2 molecules and carbon atomes describes the overall interaction in the complex rather well, and the parameters derived earlier from small hydrocarbon dopants 22 perform satisfactorily against the best electronic structure calculations performed in the present investigation for the C20 fullerene. The various sets of parameters derived from the different quantum chemical calculations for the C20 –H2 complex were transferred to C60 –H2 and, more importantly, to H2 clusters coating these fullerenes, with the aim to evaluate how sensitive are the collective behaviors on these parameters. Our investigation was here primarily concerned with the closure of the first solvation shell for which indirect measurements are available from mass spectrometry. 21 Using path-integral molecular dynamics simulations, we found that C20 and C60 are able to accomodate shells of 32 and 51–54 hydrogen molecules, respectively, a significant dependence on the potential parameters being found for the larger fullerene. Upon heating, the hydrogen shell undergoes increasingly fluxional motion with progressive melting above 5 K for C20 (H2 )32 and above 8 K for C60 (H2 )50 . For both systems, isomerization was found to set in with molecules popping away from the main shell,

21



but in the larger cluster the radial expansion makes the floater to be inserted back at the onset of melting. The nonmonotonic thermal behavior noted in C60 (H2 )50 is a manifestation of vibrational delocalization and is not expected to convey to other similar sizes, and strong finite size effects may drastically change the phenomenology of melting in neighboring complexes. In the future it would be interesting to consider much larger clusters containing multiple shells and determine whether the degrees of radial and angular fluxionality coincide or, conversely, if the outermost molecules exhibit greater mobility in a type of surface melting mechanism. Beyond the present effective pair potentials derived in the present work, more accurate anisotropic potentials could also be extracted from the same quantum chemical methods but a different functional form should then be employed, possibly accounting for the quadrupolar distribution of the hydrogen molecule. Such a more refined description could be relevant at higher densities, as in the case of condensed systems under pressure. Ionic dopants, which also interact with H2 more strongly than neutral impurities owing to polarization forces, could also benefit from a full atomistic description that no longer assumes the hydrogen molecules to freely rotate in the vicinity of a major binding site.

References (1) Liu, C.; Fan, Y. Y.; Liu, M.; Cong, H. T.; Cheng, H. M.; Dresselhaus, M. S.; Hydrogen Storage in Single-Walled Carbon Nanotubes at Room Temperature, Science 1999, 286, 11271129. (2) Cheng, H.-M.; Yang, Q.-H.; Liu, C.; Hydrogen Storage in Carbon Nanotubes, Carbon 2001, 39, 1447-1454. (3) Froudakis, G. E.; Hydrogen Storage in Nanotubes & Nanostructures, Mater. Today 2011, 14, 324-328.

22


Page 22 of 28

Page 23 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(4) Zhou, S.; Liu, X.; Zou, H.; Study of H2 Physical Adsorption in Single-Walled Carbon Nanotube Array, AIP Adv. 2013, 3, 082119. (5) Durbin, D. J.; Allan, N. L.; Malardier-Jugroot, C.; Molecular Hydrogen Storage in Fullerenes – A Dispersion-Corrected Density-Functional Theory Study, Int. J. Hydrogen Energy 2016, 41, 13116-13130. (6) Morisset, S.; Aguillon, F.; Sizun M.; Sidis, V.; Wave-packet Study of H2 Formation on a Graphite Surface Through the Langmuir-Hinshelwood Mechanism, J. Chem. Phys. 2004, 121, 6493-6501. (7) Casolo, S.; Tantardini, G. F.; Martinazzo, R.; Insights into H2 Formation in Space from Ab Initio Molecular Dynamics, Proc. Natl. Acad. Sci. USA 2013, 110, 6674-6677. (8) Pasquini, M.; Bonfanti, M.; Martinazzo, R.; Quantum Dynamical Investigation of the Isotope Effect in H2 Formation on Graphite at Cold Collision Energies, Phys. Chem. Chem. Phys. 2016, 18, 6607-6617. (9) Stan, G.; Cole, M. W.; Hydrogen Adsorption in Nanotubes, J. Low Temp. Phys. 1998, 110, 539-544. (10) Becher, M.; Haluska, M.; Hirscher, M.; Quintel, A.; Skakalova, V.; Dettlaff-Weglikovska, U.; Chen, X.; Hulman, M.; Choi, Y.; Roth, S.; et al., Hydrogen Storage in Carbon Nanotubes, Comptes Rendus Phys. 2003, 4, 1055-1062. (11) Dodziuk, H.; Modeling Complexes of H2 Molecules in Fullerenes, Chem. Phys. Lett. 2005, 410, 39-41. (12) Cross, R. J.; Does H2 Rotate Freely Inside Fullerenes? J. Phys. Chem. A 2001, 105, 69436944. (13) Riahi, S.; Pourhossein, P.; Zolfaghari, A.; Ganjali, M. R.; Jooya, H. Z.; Encapsulation of

23



Hydrogen Molecule in Fullerene (C60 ), Fullerenes, Nanotubes and Carbon Nanostructures 2009, 17, 159-170. (14) Pupysheva, O. V.; Farajian, A. A.; Yakobson, B. I.; Fullerene Nanocage Capacity for Hydrogen Storage, Nano Lett. 2008, 8, 767-774. (15) Sun, D .Y.; Liu, J. W.; Gong, X. G.;Liu, Z.-F. Empirical Potential for the Interaction between Molecular Hydrogen and Graphite, Phys. Rev. B 2007, 75, 075424. (16) Gavardi, E.; Cuppen, H. M.; Hornekær, L.; A Kinetic Monte Carlo Study of Desorption of H2 from Graphite (0001), Chem. Phys. Lett. 2009, 477, 285-289. (17) Wang, Z.; Yao, M.; Pan, S.; Jin, M.; Liu, B.; Zhang, H. A Barrierless Process from Physisorption to Chemisorption of H2 Molecules on Ligh-Element-Doped Fullerenes, J. Phys. Chem. C 2007, 111, 4473-4476. (18) Tian, C.; Wang, Z.; Jin, M.; Zhao, W.; Meng, Y.; Wang, F.; Feng, W.; Liu, H.; Ding, D.; Wu, D. Transformation Mechanism of a H2 Molecule from Physisorption to Chemisorption in Pristine and B-Doped C20 Fullerenes, Chem. Phys. Lett. 2011, 511, 393-398. (19) Arboleda Jr, N. B.; Kasai, H.; Nakanishi, H.; Scattering and Dissociative Adsorption of H2 on the Armchair and Zigzag Edges of Graphite, J. Appl. Phys. 2004, 96, 6331-6336. (20) Xu, M.; Sebastianelli, F.; Bacic, Z.; Lawler, R.; Turro, N. J.; H2 , HD, and D2 inside C60 : Coupled Translation-Rotation Eigenstates of the Endohedral Molecules from Quantum FiveDimensional Calculations, J. Chem. Phys. 2008, 128, 011101. (21) Kaiser, A.; Leidlmair, C.; Bartl, P.; Zöttl, S.; Denifl, S.; Mauracher, A.; Probst, M.; Scheier, P.; Echt, O.; Adsorption of Hydrogen on Neutral and Charged Fullerene: Experiment and Theory, J. Chem. Phys. 2013, 138, 074311. (22) Calvo, F.; Yurtsever, E.; Solvation of Carbonaceous Molecules by para-H2 and ortho-D2 Clusters. II. Fullerenes, J. Chem. Phys. 2016, 145, 084304. 24


Page 24 of 28

Page 25 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(23) Calvo, F.; Yurtsever, E.; Solvation of Carbonaceous Molecules by para-H2 and ortho-D2 Clusters. I. Polycyclic Aromatic Hydrocarbons, J. Chem. Phys. 2016, 144, 224302. (24) Buck, U.; Huisken, F.; Kohlhase, A.; Otten, D.; Schaefer, J.; State Resolved Rotational Excitation in D2 +H2 Collisions, J. Chem. Phys. 1983, 78, 4439-4450. (25) Silvera, I. F.; Goldman, V. V.; The Isotropic Intermolecular Potential for H2 and D2 in the Solid and Gas Phases, J. Chem. Phys. 1978, 69, 4209-4213. (26) Schelkacheva, T. I.; Tareyeva, E. E.; Orientational Phase Transition in Solid C60 , Phys. Rev. B 2000, 61, 3143-3146. (27) Pérez, A.; Tuckerman, M. E.; Müser, M. H. A.; A Comparative Study of the Centroid and Ring-Polymer Molecular Dynamics Methods for Approximating Quantum Time Correlation Functions from Path Integrals, J. Chem. Phys. 2009, 130, 184105. (28) Witt, A.; Ivanov, S. D.; Forbert, H.; Marx, D.; On the Applicability of Centroid and Ring Polymer Path Integral Molecular Dynamics for Vibrational Spectroscopy, J. Chem. Phys. 2009, 130, 194510. (29) Mezzacapo, F.; Boninsegni, M.; On the Possible “Supersolid” Character of Parahydrogen Clusters, J. Phys. Chem. A 2011, 115, 6831-6837. (30) Hobza, P.; Muller-Dethlefs, K.; Non-covalent Interactions: Theory and Experiment, RSC Publishing, London (2009). (31) Grimme, S.; Improved Second-Order Møller-Plesset Perturbation Theory by Separate Scaling of Parallel- and Antiparallel-spin Pair Correlation Energies, J. Chem. Phys. 2003, 118, 90959102. (32) Jung, Y.; Lochan, R. C.; Dutoi, A. D.; Head-Gordon, M.; Scaled Opposite-Spin Second Order Møller-Plesset Correlation Energy: An Economical Electronic Structure Method, J. Chem. Phys. 2004, 121, 9793-9802. 25



(33) Distasio Jr, R. A.; Head-Gordon, M. Optimized Spin-Component Scaled Second-Order Møller-Plesset Perturbation Theory for Intermolecular Interaction Energies, Mol. Phys. 2007, 105, 1073-1083. (34) Jeziorski, B.; Moszynski, R.; Szalewicz, K.; Perturbation Theory Approach to Intermolecular Potential Energy Surfaces of van der Waals Complexes, Chem. Rev. 1994, 94, 1887-1930. (35) Szalewicz, K.; Patkowski, K.; Jeziorski, B. Intermolecular Interactions via Perturbation Theory: From Diatoms to Biomolecules, Struct. Bonding. 2005, 116, 43-117. (36) Jansen, G.; Heßelmann, A.; Comment on “Using Kohn-Sham Orbitals in Symmetry-Adapted Perturbation Theory to Investigate Intermolecular Interactions”, J. Phys. Chem. A 2001, 105, 11156-11157. (37) Heßelmann, A.; Jansen, G.; First-Order Intermolecular Interaction Energies from the KohnSham Orbitals, Chem. Phys. Lett. 2002, 357, 464-470. (38) Heßelmann, A.; Jansen, G.; Intermolecular Induction and Exchange-Induction Energies from Coupled-Perturbed Kohn-Sham Density Functional Theory, Chem. Phys. Lett. 2002, 362, 319-325. (39) Heßelmann, A.; Jansen, G.; Intermolecular Dispersion Energies from Time-Dependent Density Functional Theory, Chem. Phys. Lett. 2003, 367, 778-784. (40) Heßelmann, A.; Jansen, G.; Schütz, M.; Density-Functional Theory-Symmetry-Adapted Intermolecular Perturbation Theory with Density Fitting: A new Efficient Method to Study Intermolecular Interaction Energies, J. Chem. Phys. 2005, 122, 014103. (41) Werner, H. J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; Celani, P.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; et al., Molpro version 2009.1, a package of ab initio programs, 2009, see http://www.molpro.net, (Accessed may 8, 2017).

26


Page 26 of 28

Page 27 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(42) Weigend, F.; A Fully Direct RI-HF Algorithm: Implementation, Optimised Auxiliary Basis Sets, Demonstration of Accuracy and Efficiency, Phys. Chem. Chem. Phys. 2002, 4, 42854291. (43) Weigend, F.; Köhn, A.; Hättig, C. Efficient Use of the Correlation Consistent Basis Sets in Resolution of the Identity MP2 Calculations, J. Chem. Phys. 2002, 116, 3175-3183. (44) Boys, S. F.; Bernardi, F. The Calculation of Small Molecular Interactions by the Differences of Separate Total Energies. Some Procedures with Reduced Errors, Mol. Phys. 1970, 19, 553-566. (45) Bak, K. L.; Jørgensen, P.; Olsen, J.; Helgaker, T.; Klopper, W.; Accuracy of Atomization Energies and Reaction Enthalpies in Standard and Extrapolated Electronic Wave Function/Basis Set Calculations, J. Chem. Phys. 2000, 112, 9229-9242. (46) Tekin, A.; Jansen, G.; How Accurate is the Density Functional Theory Combined with Symmetry-Adated Perturbation Theory Approach for CH-π and π–π Interactions? A Comparison to Supermolecular Calculations for the Acetylene-Benzene Dimer, Phys. Chem. Chem. Phys. 2007, 9, 1680-1687. (47) Zhou, Y.; Karplus, M.; Ball, K. D.; Berry, R. S.; The distance fluctuation criterion for melting: Comparison of square-well and Morse potential models for clusters and homopolymers, J. Chem. Phys. 2002, 116, 2323-2329. (48) Calvo, F. Size-Induced Melting and Reentrant Freezing in Fullerene-Doped Helium Clusters, Phys. Rev. B 2012, 85, 060502(R).

27



Graphical TOC Entry 50.0

Shell size

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 28

49.8 49.6 49.4 49.2 49.0 0

4

8

Temperature (K)

28


12

by Hydrogen Clusters at

Recommend Documents