The ground state and electronic structure of Gd@C82: a systematic theoretical investigation of first principle density functionals.

The ground state and electronic structure of Gd@C82: A systematic theoretical investigation of first principle density functionals Xing Dai, Yang Gao, Minsi Xin, Zhigang Wang, and Ruhong Zhou Citation: The Journal of Chemical Physics 141, 244306 (2014); doi: 10.1063/1.4904389 View online: http://dx.doi.org/10.1063/1.4904389 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/141/24?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Systematic theoretical investigation of the zero-field splitting in Gd(III) complexes: Wave function and density functional approaches J. Chem. Phys. 142, 034304 (2015); 10.1063/1.4905559 An accurate first principles study of the geometric and electronic structure of B 2 , B 2 − , B 3 , B 3 − , and B 3 H : Ground and excited states J. Chem. Phys. 132, 164307 (2010); 10.1063/1.3389133 An ab initio investigation on the endohedral metallofullerene Gd 3 N – C 80 J. Appl. Phys. 101, 09E105 (2007); 10.1063/1.2711420 Electronic ground states of the V 2 O 4 +/0/− species from multireference correlation and density functional studies J. Chem. Phys. 120, 4207 (2004); 10.1063/1.1643891 A density-functional study of the structures and electronic properties of C 59 Ni and C 60 Ni clusters J. Chem. Phys. 114, 9371 (2001); 10.1063/1.1353583

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 137.30.242.61 On: Fri, 08 May 2015 08:51:27

THE JOURNAL OF CHEMICAL PHYSICS 141, 244306 (2014)

The ground state and electronic structure of Gd@C82: A systematic theoretical investigation of first principle density functionals Xing Dai,1 Yang Gao,1 Minsi Xin,1,2 Zhigang Wang,1,3,a) and Ruhong Zhou1,4,5,b)

1

Institute of Atomic and Molecular Physics, Jilin University, Changchun 130012, People’s Republic of China School of Science, Changchun University of Science and Technology, Changchun 130022, People’s Republic of China 3 State Key Laboratory of Theoretical and Computational Chemistry, Jilin University, Changchun 130023, People’s Republic of China 4 Computational Biology Center, IBM Thomas J. Watson Research Center, Yorktown Heights, New York 10598, USA 5 Department of Chemistry, Columbia University, New York, New York 10027, USA 2

(Received 16 September 2014; accepted 4 December 2014; published online 24 December 2014) As a representative lanthanide endohedral metallofullerene, Gd@C82 has attracted a widespread attention among theorists and experimentalists ever since its first synthesis. Through comprehensive comparisons and discussions, as well as references to the latest high precision experiments, we evaluated the performance of different computational methods. Our results showed that the appropriate choice of the exchange-correlation functionals is the decisive factor to accurately predict both geometric and electronic structures for Gd@C82. The electronic structure of the ground state and energy gap between the septet ground state and the nonet low-lying state obtained from pure density functional methods, such as PBE and PW91, are in good agreement with current experiment. Unlike pure functionals, the popularly used hybrid functionals in previous studies, such as B3LYP, could infer the qualitative correct ground state only when small basis set for C atoms is employed. Furthermore, we also highlighted that other geometric structures of Gd@C82 with the Gd staying at different positions are either not stable or with higher energies. This work should provide some useful references for various theoretical methodologies in further density functional studies on Gd@C82 and its derivatives in the future. C 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4904389]

I. INTRODUCTION

Endohedral metallofullerenes (EMFs) are so-called “core-shell” molecules formed by encapsulating metal atoms/ clusters into various hollow carbon fullerene cages. Due to the embedded metals, significant charge transfer often occurs between the inner species and the cage which leads to the novel structures and properties for EMFs.1,2 From the theoretical point of view, it is unrealistic to perform high-precision ab initio calculations for EMFs systems based on the current available computing power, especially for lanthanide or actinide EMFs due to the large amount of electrons contained. Density functional theory (DFT) with relativistic effective core potential (RECP) approach is an effective method to describe EMFs systems and is widely used in relevant studies.3 A significant advantage of DFT calculations is that it accounts for the electron correlation.4–6 The RECP approach can effectively consider the scalar relativistic effects and greatly save computational time7,8 which is crucial for reasonably describing the electronic structures of EMFs. Previous studies have reported that the spin-orbit coupling (SOC) effects do not significantly affect the equilibrium structure and the ground state.9–12 Meanwhile, this effect is generally quenched in large molecules13 such as EMF.14 Therefore, generally the spin-orbit a)Electronic mail: [email protected] b)Electronic mail: [email protected]

0021-9606/2014/141(24)/244306/12/$30.00

coupling effect is omitted in theoretical studies of EMFs. DFT with RECP approach has become the most popular method in computational EMFs studies with calculation results in good agreement with the experimental observations in general. However, due to the lack of comprehensive DFT assessments on EMFs systems, theoreticians often choose DFT methods empirically to deal with their EMFs systems in hand, including the choices of exchange-correlation (XC) functionals and basis sets. Consequently, (slightly) different computational results are obtained which often raise arguments. For example, Gd@C82 is one of the most controversial representatives. Because of the important potential applications of Gd@C82 and its derivatives in the biomedical field,15–18 they have been widely investigated in the past two decades. Nevertheless, both the stable geometry and ground state of Gd@C82 still remain controversial in these relevant works. In the following paragraphs, we summarize and review these studies from both experimental and theoretical aspects. Funasaka et al. isolated a macroscopic quantity of the EMF Gd@C82 in 1995. The magnetization data obtained employing a Quantum Design SQUID magnetometer showed a signature of paramagnetic behavior and inferred the gadolinium as a Gd3+ free ion in fullerene.19,20 Suzuki et al. performed the cyclic voltammogram (CV) studies and suggested that the oxidation states of the Gd in C82 fullerene are 3+, and the 4f electrons of Gd do not play any important role in fullerenolanthanide chemistry.21 The valence state expression

141, 244306-1

© 2014 AIP Publishing LLC


244306-2

Dai et al.

of Gd3+@C823− has been confirmed in the transmission electron energy-loss spectroscopy (EELS) experiment.22 Subsequently, Huang et al. studied the magnetic properties of Gd@ C82 and revealed a septet ground spin state, where the energy gap between the ground state (septet) and the higher spin state (nonet) is more than 25 meV.23,24 However, this value was considered as temperature-independent, which needs further hightemperature susceptibility measurement to clarify.25 Then, Furukawa et al. performed the first quantitative estimate of the magnetic interaction between the metal ion and radical spin on the cage by a high-field electron spin resonance (ESR) spectrometer and reported that the ground spin state of Gd@C82 is resulted from the antiferromagnetically intramolecular exchange coupling between SGd = 7/2 (the 7 net spin electrons on Gd) and SR = 1/2 (the 1 net spin electrons on the cage). Also, the further spectrum of the Gd@C82 monomer was assigned to the mixed spin states of the ground S = 3 (septet) and the thermally excited S = 4 (nonet) states with an energy gap of 14.4 cm−1 (1.79 meV)26 which provides an important reference for the subsequent studies. An extended x-ray absorption fine structure (EXAFS) study27 has shown that the nearest neighbor distance from the encaged Gd atom to carbon atom is 2.49(3) Å, and the results can be well fitted with a Gd positioned above a carbon hexagon, indicating Gd@C82 should have a normal M@C82 geometric structure similar to Sc@C82 and [email protected],29 The cage symmetry of Gd@C82 has been speculated to be C2v by UV/Vis/NIR spectra.30 However, Nishibori performed an experiment of synchrotron radiation powder structure analysis by using the maximum entropy method (MEM) and revealed that Gd@C82 has an anomalous geometric structure, in which the Gd atom is located in the vicinity of the C–C bond on the C2 axis in the C82 (C2v ) cage, which is the opposite end of the six-membered ring where Sc and La atoms in Sc@C82 and La@C82 are known to be located.31 With the continuous advancement of experiment technology in recent years, both the X-ray absorption near-edge structure (XANES) study32 and the single crystal X-ray structure of the Carbene Adduct experiment33 showed that Gd@C82 should be a normal M@C82 structure with the Gd atom near the sixmembered ring on the C2 axis. The latest single-crystal X-ray diffraction analyses of the cocrystals of Gd@C82 with nickelII octaethylporphyrin [NiII (OEP)] further confirmed the normal Gd@C82 geometry.34 On the theoretical side, different calculation methods tend to generate different stable geometries and electronic states of Gd@C82. In the earliest theoretical study, with the fixed C82 cage, B3LYP functional predicted a nonet ground state arising from the ferromagnetic coupling between the unpaired electrons of the Gd-4f and the cage.35 Coincidentally, another B3LYP full-optimization study also predicted a nonet ground state and revealed an anomalous geometric structure of Gd@C82 with the Gd lying on the C2 axis and close to the C–C bond which agreed with the MEM experimental observation.36 For some period of time, the stable Gd@C82 was considered to have a nonet ground state but with a more anomalous geometry where the Gd atom deviates from the C2 axis and adsorbs another C–C bond.37,38 Consequently, more detailed DFT-B3LYP studies considering various Gd@C82 isomers showed that, in fact, the so-called anomalous structure

J. Chem. Phys. 141, 244306 (2014)

is not stable, and the most stable Gd@C82 should adopt a normal M@C82 atomic configuration. Although this report agreed with the experimental observed septet ground state, the theoretical calculations still revealed a nonet ground state which is only 0.01 kcal/mol (corresponding to 0.43 meV) lower than septet.32,39 Additionally, a theoretical study of M@C82 containing other metals also indirectly supports that Gd should locate on the C2-hexagon in C82 cage.40 The latest study via generalized gradient approximation (GGA) plus Hubbard U model (GGA+U) of DFT explained the origin of the antiferromagnetic coupling between endohedral Gd and the free spin on the carbon cage and identified a 10 meV septetnonet gap.25 In Ref. 34, the authors also performed M06-2X calculations and confirmed that the most favorable site of the inner Gd is in agreement with their experimental result, but they did not discuss the electronic structure of the ground state. Although it appears that the controversies surrounding the geometry and electronic state of Gd@C82 have been resolved, the effects of different theoretical methods will stay for a long time. This is true especially for the identification of the ground state since there is no reliable method for determining this due to the lack of reliable and systematic evaluation for DFT. Therefore, to help the exploration of the electronic structure of unknown EMFs, it is necessary to determine more reliable theoretical methods. In this work, we performed a systematic study for Gd@C82 by considering the two aspects of XC functionals and basis set aiming at establishing a regularity that clearly reflects their effects on EMFs under the DFT framework. Gd@C82 is chosen here mainly because it has been widely explored, with some methods consistent with the latest experimental data, which makes it the best candidate for benchmarking. The advantages and disadvantages of different calculation methods can be observed from the comprehensive results, and finally we summarized the relatively reliable and suitable theoretical approaches for this EMF. These results are vital for the subsequent studies of Gd@C82, as well as its derivatives. II. COMPUTATIONAL DETAILS

We performed full geometric optimization for Gd@C82 using ten standard DFT methods in Gaussian 09 program.41 To take into account the scalar relativistic effects, the effective core potential (ECP) triple split basis set (CEP-121G)42 was applied for the Gd atom which is used in almost all previous theoretical studies of Gd@C82. For the CEP-121G basis set, the Gd’s all possible active electrons (i.e., 4f, 5s, 5p, 5d, and 6s electrons, corresponding to 18 valence electrons) are allowed to relax in the electronic structure calculations, and the relativistic effect is implicitly described in the corresponding pseudo potential part.42 Therefore, it is a good candidate for both computational cost and accuracy in lanthanide studies.32,36,43 Meanwhile, a more accurate treatment of 28 core-electron pseudo potential44 developed by Dolg et al. with its corresponding [14s13p10d8f6g]/(6s6p5d4f3g) valence basis set45 was also used for Gd, and correlated calculated results was presented in the supplementary material (Ref. 90) part 2. The SOC was not considered because it does not significantly affect the equilibrium structure and ground


244306-3

Dai et al.

state.9–12 Previous report has pointed out that SOC effects are diluted by the presence of the fullerene cage, and there are no large effects on the ground-state molecular properties.46 In fact, more sophisticated treatments of relativistic effects (both the scalar part and SOC) are generally highly resourcedemanding and time-consuming making them unaffordable for such large system like Gd@C82. Based on these considerations, we employed the common DFT/RECP method to ensure that the electron correlations and the most significant relativistic effects are covered in the present study. In part 3 of supplementary material (Ref. 90), we also discussed the results calculated at DKH-DFT (all-electron scalar relativistic DFT method utilizing the Douglas-Kroll-Hess Hamiltonian) and ZORA-DFT (scalar relativistic DFT method using Zero Order Regular Approximation) methods. Our previous works showed that the nature of the basis sets of carbon atoms should impact the electronic structures of spin-polarization carbon materials.47–52 In this study, we were also concerned about the impacts on the electronic structure of Gd@C82 from the basis set. For the sake of comparability, universality, and computing capacity, the 3-21G,53 4-31G,54 4-31G*,54 6-31G,55 6-31G*,55,56 6-31+G*,55–57 6-311G,58 6-311G*,58 and 6-311+G*57,58 basis sets were chosen for the C atoms. These choices will help us systematically determine the influences from the basis sets, including their size, the introduction of polarization or diffuse functions, and the split forms (double or triplet split). As an additional validation, the Ahlrichs59,60 (def2-SVP and def2-TZV) and Dunning (ccpvdz,61 cc-pvtz62) types of basis sets were also employed for C atoms. Both the local spin density approximation (LSDA)63–66 and GGA were employed. For the GGA method, the PBE,67,68 PW91,69–73 BLYP,74–76 BP86,74,77 BPW91,69–74 B3LYP,75,76,78 B3PW91,69–73,78 and PBE079 functionals were chosen. We also considered the dispersion effects by utilizing the third generation dispersion-corrected DFT-D380 including the PBE-D3 and B3LYP-D3 methods. The hybrid meta XC functional M062X81 was also explored. Through our tests, the quintet and 11-et states of Gd@C82 always have much higher total energies (more than 1 eV) than the septet and nonet states. This conclusion was also presented in previous reports.32 Therefore in this work, we only performed full geometric optimization calculations for the septet and nonet states of Gd@C82 and focused on their total energies. For the PBE, PW91, B3LYP, PBE0, and B3PW91 functionals, the frequency calculations were also performed after the geometric optimization in order to ensure the structures we obtained are real minima on the potential energy surface and consider the zero point energy (ZPE) corrections. III. RESULTS

Through the descriptions above, some primary facts have to be realized before the discussion. In the C82 (C2v ) fullerene, the C2 axis crosses the center of a hexagon and a C–C bond. In fact, the neutral C82 (C2v ) is unstable due to its open-shell electronic structure which cannot be obtained as an empty cage.82 A possibility of obtaining C82 (C2v ) is from endohedral metallofullerenes,82 such as Gd@C82. After years of debates,

J. Chem. Phys. 141, 244306 (2014)

FIG. 1. The diagram of Gd@C82 (C 2v ) structure. The C 2 axis crosses the centers of a hexagon (top, yellow) and a C–C bond (bottom, blue). The Gd atom (green) lies on the C 2 axis and near the top hexagon. This type of atomic configuration is derived from the latest experimental observation and is commonly known as “normal M@C82.” This geometry of Gd@C82 is also in accordance with our predicted structure with global energy minimum.

the most stable atomic configuration of Gd@C82 is finally determined by experiments as a normal M@C82 structure with the Gd atom staying at the C2 axis as described in Fig. 1. According to the existing experimental results, the nearest GdC distance is about 2.49(3) Å (Ref. 27), and the septet-nonet gap is about 1.79 meV.26 Because the energy gap is so small, the septet and nonet states were often confused in previous DFT theoretical studies. The experimentally determined geometric parameters and the septet-nonet gap can be used to test the quality of XC functionals. The discussion of our results is divided into three parts. First, we performed full geometric optimization calculations for the experimentally confirmed Gd@C82 structure (Fig. 1) using different DFT methods to evaluate the reliability of XC functionals and basis sets of C. Second, we explain the ground state electronic structure of Gd@C82 based on the reliable methods evaluated. Third, we further check other geometries with Gd located at other sites in C82 cage which have been widely disputed before. A. Evaluations on exchange-correlation functionals and the basis sets of carbon atoms 1. Pure functionals PBE, PW91, BLYP, BP86, and BPW91

We first compare the results calculated at two pure GGA functionals, namely, PBE and PW91. Table I indicates that both the PBE and PW91 functionals correctly distinguish the ground state of Gd@C82 by figuring out that the septet state has the lowest total energy. The septet-nonet gaps are very close (at meV magnitude) from these two methods and are in good agreement with the ESR experiment.26 A series of rules about the influences of the carbon’s basis sets on the geometric parameters and septet-nonet energy gap can be summarized as follows: (I) If neither polarization nor diffuse function is introduced, the energy gap will reduce when only the scale of the carbon’s basis sets is enlarged (compare the results of 321G, 4-31G, and 6-31G). When the triplet spit 6-311G basis set is used, the energy gap is slightly larger than the double split 631G result. (II) For the double split basis sets, the introduction


244306-4

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

TABLE I. The septet-nonet gap and Gd-C distances of Gd@C82 calculated at pure PBE and PW91 levels. ∆E = E(nonet) − E(septet). (Note: Frequency calculations with large basis set contained are extremely timeconsuming, so it is not performed for some cases.) The boldface values indicate that they are closed to experimental value, with or without ZPE correction. Gd-C (Å) Expt. = 2.49(3)

Method EXC PBE

PW91

Basis set

∆E (meV) (+ZPE) Expt. = 1.79

Septet

Nonet

3.23 (2.59) 2.13 (1.55) 2.32 (1.88) 1.60 (1.03) 1.87 (1.44) 1.28 (1.09) 1.90 (0.68) 1.90 (1.61) 1.85 1.44 (1.09) 2.56 2.24 (1.85) 2.66 3.14 (2.50) 2.06 (1.61) 2.25 (1.80) 1.55 (1.12) 1.80 (1.36) 1.28 (1.11) 1.84 (0.44) 1.83 (1.44) 1.79 1.42 (1.09) 2.48 2.14 (1.74) 2.61

2.48, 2.49 2.48, 2.48 2.46, 2.47 2.48, 2.49 2.46, 2.47 2.46, 2.47 2.48, 2.49 2.46, 2.46 2.47, 2.47 2.46, 2.46 2.45, 2.45 2.47, 2.47 2.46, 2.46 2.48, 2.49 2.48, 2.48 2.46, 2.47 2.48, 2.49 2.46, 2.47 2.46, 2.47 2.48, 2.49 2.46, 2.46 2.47, 2.47 2.46, 2.46 2.45, 2.45 2.47, 2.47 2.46, 2.46

2.48, 2.49 2.48, 2.48 2.46, 2.47 2.48, 2.49 2.46, 2.47 2.46, 2.47 2.48, 2.49 2.46, 2.46 2.47, 2.47 2.46, 2.46 2.45, 2.45 2.47, 2.47 2.46, 2.46 2.48, 2.49 2.48, 2.48 2.46, 2.47 2.48, 2.49 2.46, 2.47 2.46, 2.47 2.48, 2.49 2.46, 2.46 2.47, 2.47 2.46, 2.46 2.45, 2.45 2.47, 2.47 2.46, 2.46

CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼def2-TZVP CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼def2-TZVP

of the d polarization functions generally increases the energy gap (about 0.2–0.3 meV, compare the results of 4-31G and 4-31G*, 6-31G and 6-31G*). And for the triplet split basis sets, the energy gap is almost not affected by the polarization functions (compare the results of 6-311G and 6-311G*). (III) For double and triplet split basis, the energy gap will be reduced by about 0.5 meV and 0.05 meV, respectively, when the diffuse functions are introduced (compare the results of 6-31G* and 6-31+G*, 6-311G* and 6-311+G*). (IV) For double split basis, the ZPE corrections will generally reduce the energy gap by about 0.4–0.6 meV. This is because the ZPE correction (positive number) for the septet is often larger than that for the nonet, and consequently the total energies of the septet and nonet become closer. However, when considering the ZPE corrections, the result of triplet split 6-311G basis set (without polarization function) seriously underestimates the energy gap. For the result of triplet split 6-311G* basis set (with d polarization function), the energy gap is reduced by 0.3–0.4 meV, which can be compared with the 6-31G* results. (V) In terms of the geometric parameters, there is almost no effect when only the size of the carbon’s basis sets is changed (compare Gd-C distances of 3-31G, 4-31G, 6-31G, and 6-311G). The calculated Gd-C distances are about 2.48–2.49 Å which is in excellent agreement with the experimental observation (2.49(3) Å). The introduction of polarization functions will bring about

0.01–0.03 Å, reduction for the Gd-C distance, and the diffuse functions will result in an increase of about 0.01 Å. The results calculated at BLYP, BP86, and BPW91 (all of them employ the B exchange functional) are listed in Table II. For BLYP, an obvious phenomenon is that it always significantly overestimates the energy gap when polarization function is introduced except for 6-31G*. That is, because BLYP calculates the stable geometry with longer Gd-C distance (2.62-2.66 Å, Table II) and overestimates the stability of the septet when the carbon’s basis set contains the d polarization function. This is due to BLYP determined an incorrect septet electronic arrangement for the geometry with longer Gd-C distance (here, the net spins of Gd and cage are about 6 (spin up) and 0 (spin down), respectively). Therefore, BLYP presents uncertainties in determining the energy gap. For BP86 and BPW91, all calculations show a small energy gap, while the values of BP86 can be smaller than 1 meV when C’s basis sets become larger. Here, the effects from carbon’s basis set are similar to that of PBE and PW91. Summarizing Table I and Table II, all results of these pure XC functionals could qualitatively support the ground state of Gd@C82 being septet. The calculated septet-nonet gap of PBE and PW91 agrees well with the ESR experimental observation and does not significantly depend on the C’s basis set. In addition, different C’s basis sets will not result in miscalculation on


244306-5

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

TABLE II. The septet-nonet gap and Gd-C distances of Gd@C82 calculated at pure BLYP, BP86, and BPW91 levels. ∆E = E(nonet) − E(septet). The boldface values indicate that they are closed to experimental value, with or without ZPE correction. Gd-C (Å) Expt. = 2.49(3)

Method EXC BLYP

BP86

BPW91

Basis set

∆E (meV) Expt. = 1.79

Septet

Nonet

CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP

2.55 1.36 861.65 0.79 1.36 903.36 1.17 898.71 939.09 863.85 1.73 1.38 2.32 1.18 1.46 0.68 0.98 0.49 0.93 0.97 1.02 0.62 1.52 1.34 2.84 1.76 1.95 1.24 1.49 1.14 1.53 1.49 1.62 1.11 2.13 1.89

2.50, 2.50 2.50, 2.50 2.64, 2.62 2.50, 2.51 2.49, 2.49 2.65, 2.63 2.51, 2.51 2.65, 2.63 2.66, 2.64 2.63, 2.62 2.48, 2.48 2.49, 2.49 2.49, 2.49 2.48, 2.49 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.47 2.48, 2.49 2.47, 2.47 2.47, 2.47 2.47, 2.47 2.46, 2.46 2.47, 2.47 2.48, 2.49 2.48, 2.49 2.46, 2.47 2.48, 2.49 2.47, 2.47 2.47, 2.47 2.48, 2.49 2.46, 2.47 2.47, 2.47 2.46, 2.47 2.45, 2.46 2.47, 2.47

2.50, 2.50 2.50, 2.50 2.49, 2.49 2.50, 2.51 2.49, 2.49 2.49, 2.49 2.51, 2.51 2.49, 2.49 2.49, 2.49 2.48, 2.49 2.48, 2.48 2.49, 2.49 2.49, 2.49 2.48, 2.49 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.47 2.47, 2.47 2.46, 2.46 2.47, 2.47 2.48, 2.49 2.48, 2.49 2.46, 2.47 2.48, 2.49 2.47, 2.47 2.47, 2.47 2.48, 2.49 2.46, 2.47 2.47, 2.47 2.46, 2.47 2.45, 2.46 2.47, 2.47

the ground state when using the pure functional DFT methods. From these results, the double-split 4-31G* and 6-31G* basis sets are relatively good choices for the C atoms from both accuracy and computational efficiency point of views. In fact, even the small 3-21G basis set is chosen, the obtained septetnonet gap and Gd-C distance are still acceptable. Our previous study has shown that the reliability of the electronic state is not simply determined by the size of the basis sets.49 There is no doubt that the choice of small basis sets could significantly improve the computational efficiency, therefore it is worth considering smaller C’s basis sets when performing calculations without high-precision requirement. 2. Hybrid functionals B3LYP, PBE0, and B3PW91

Hybrid exchange functionals, such as the Becke three parameter hybrid functional (B3), are parameterized using a

test set which contains no transition metals, actinides or lanthanide.75,76,78 Previous study showed that, in transition metal dimer complex, B3LYP performs poorly for both structural and energy gap (between ground state and low-lying excited state) predictions.83 Therefore, such functionals might not be expected to perform well on Gd@C82 which has a fairly small septet-nonet gap. In fact, the electronic structure of U@C82 predicted by B3LYP is in agreement with the experiment,84 but this molecule does not involve the issue of small energy gap between the respective electronic states. For Gd@C82, B3LYP functional was used in almost all previous theoretical studies and did not give the correct ground state. Therefore, such small septet-nonet gap may be sensitive to the hybrid functionals. In this section, we summarized the energy and geometric parameters of Gd@C82 calculated with B3LYP, as well as PBE0 and B3PW91 functionals, in Table III, and systematically evaluated the feasibility of these methods.


244306-6

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

TABLE III. The septet-nonet gap and Gd-C distances of Gd@C82 calculated at hybrid B3LYP, PBE0, and B3PW91 levels. ∆E = E(nonet) − E(septet). Gd-C (Å) Expt. = 2.49(3)

Method EXC B3LYP

PBE0

B3PW91

Basis set

∆E (meV) (+ZPE) Expt. = 1.79

Septet

Nonet

1.17 (1.03) −0.27(−0.68) −0.09(−0.44) −0.93(−1.36) −0.65(−0.93) −1.55 −0.74 −0.76 −0.95 −1.08 −0.26 −0.34 −0.32 1.78 (1.12) 0.33 (−0.19) 0.28 (−0.08) −0.31(−0.79) −0.28(−0.60) 2.05 (1.42) 0.67 (0.11) 0.67 (0.24) 0.04 (−0.49) 0.11 (−0.22)

2.49, 2.49 2.48, 2.49 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.48 2.47, 2.47 2.46, 2.46 2.47, 2.48 2.47, 2.47 2.47, 2.48 2.47, 2.48 2.45, 2.46 2.47, 2.48 2.45, 2.46 2.47, 2.48 2.47, 2.48 2.45, 2.46 2.47, 2.48 2.45, 2.46

2.49, 2.49 2.48, 2.49 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.47 2.49, 2.49 2.47, 2.47 2.47, 2.48 2.47, 2.47 2.46, 2.46 2.47, 2.48 2.47, 2.47 2.47, 2.48 2.47, 2.47 2.45, 2.46 2.47, 2.48 2.45, 2.46 2.47, 2.48 2.47, 2.48 2.45, 2.46 2.47, 2.48 2.45, 2.46

CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼def2-TZVP CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G*

Obviously, the results of the septet-nonet gap calculated at B3LYP are strongly dependent on the choice of C’s basis set. Interestingly, only the 3-21G basis set could infer the correct ground state as septet for the B3LYP functional. Once the C’s basis set is enlarged, B3LYP always gives a wrong conclusion. This phenomenon may well explain the reason why a nonet ground state was identified by previous B3LYP studies.32,35–37,39 The PBE0 results demonstrate that enlarging the C’s basis set would result in closer energy gaps, even with a wrong energy order. B3PW91 could obtain the qualitatively right energy order using these basis sets selected in Table III, but the energy order would also reverse after the ZPE correction. From these discussions, it can be inferred that hybrid functionals are unreliable for calculating the energy gap between different electronic states since they often poorly predict the energy gap and even lead to a wrong energy order. Therefore, extreme care must be taken when using hybrid DFT and C’s basis sets for the electronic state studies of similar EMFs in the future, especially for the possibility of small energy gap. 3. LSDA

Despite the fact that DFT-GGA method is currently widely accepted in EMFs studies while LSDA is not popularly used, the quality of calculating EMFs of LSDA still needs to be fully understood. Table IV shows that LSDA always gets a qualitatively correct septet-nonet energy order without fully depending on the C’s basis sets. However, the energy gap is larger than the experimental result (about twice larger) as well

as PBE and PW91 results. Also, the Gd-C distance is slightly shorter than the experimental value. 4. M06-2X and DFT-D3

Usually, the obvious metal-cage (M-C) ionic interaction (electrostatic interaction) can be attributed to the electrons donated by the inner metals to the carbon cage in EMFs.85 Meanwhile, the d or f atomic orbitals of the metals could hybridize with the π orbitals of the cage which also contribute to the covalent feature.85 Since all these ionic M-C, covalent M-C, and C–C interactions are quite strong, previous EMFs studies generally did not consider the weak interaction such as the dispersion effect. However, for the special Gd@C82 which has a small energy gap between the ground state and lowlying exited state, the dispersion correction may also affect the septet-nonet gap or order. Thus, we also performed M06-2X and DFT-D3 calculations for Gd@C82. The M06-2X functional contains dispersion effect in its Hamiltonian and takes part in the self consistent field (SCF) iteration calculations. This functional was recommended for applications involving main-group thermochemistry, kinetics, and noncovalent interactions.81 In the latest Gd@C82 study,34 the authors theoretically confirmed the most stable position of Gd by using M06-2X functional, but they did not discussed the electronic state or energy gap. From Table V, no matter what basis set is chosen for C atoms, M06-2X could always predict the septet ground state but with a large septet-nonet gap (about 300-400 and 800-1000 meV). This is because both the septet and nonet electronic arrangements obtained from M06-


244306-7

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

TABLE IV. The septet-nonet gap and Gd-C distances of Gd@C82 calculated at LSDA. ∆E = E(nonet) − E(septet). Gd-C distance (Å) Expt. = 2.49(3)

Method EXC LSDA

Basis set

∆E (meV) Expt. = 1.79

Septet

Nonet

CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP

5.35 4.27 4.28 3.61 3.70 2.96 4.27 3.89 3.94 3.33 4.72 4.33

2.45, 2.45 2.44, 2.45 2.43, 2.43 2.44, 2.45 2.43, 2.43 2.43, 2.43 2.44, 2.46 2.43, 2.43 2.43, 2.43 2.42, 2.43 2.42, 2.42 2.43, 2.43

2.45, 2.45 2.44, 2.45 2.43, 2.43 2.44, 2.45 2.43, 2.43 2.43, 2.43 2.44, 2.46 2.43, 2.43 2.43, 2.43 2.42, 2.43 2.42, 2.42 2.43, 2.43

2X are incorrect (here, the net spins of Gd is about 6 (spin up), the net spins of cage is about 0 (spin down) or 2 (spin up)). The Gd-C distances calculated at M06-2X are also longer than the experimental data and other GGA results (about 0.2 Å longer). In fact, the authors of Ref. 81 have pointed out that the M06-2X functional is a high-nonlocality functional with double amount of nonlocal exchange (2X), and it is parameterized only for nonmetals excluding the lanthanides. Combined with the results in Table V, the M06-2X functional is likely unsuitable for Gd@C82. Dispersion-corrected DFT-D is another method to consider the dispersion effects developed by Grimme.80 In this method, dispersion correction is treated as an add-on to the self-consistent Kohn-Sham energy obtained from the chosen standard DFT. The DFT-D method is suitable for the first 194 elements in the periodic table. Here, the latest generation of DFT-D3 method was employed to consider the dispersion effect for Gd@C82, and the results are listed in Table VI. Comparing Table I and Table VI, the geometric parameters from DFT-D3 are very close to that from the pure func-

tionals. Meanwhile, the geometries of the septet and nonet are also similar. According to the expression of DFT-D3,80 this dispersion-corrected method does not modify the Kohn-Sham potential in essence and the dispersion energy only depends on the geometric structure. Therefore, the dispersion energies of the septet and nonet are as close as the ones seen in Table VI. For PBE-D3, the dispersion energy is about 2.96 eV (accounted for about 0.3% of the total energy). The correction (dispersion) energy of nonet is 0.15 meV larger than septet. As a result, the septet-nonet gap from PBE-D3 reduces by about 0.15 meV compared to the standard PBE (Table I). Comparing Table III and VI, however, the dispersion correction of nonet is also larger than septet at B3LYP-D3, so the introduction of the dispersion correction could not change the qualitative conclusion on the septet-nonet energy order compared to the standard B3LYP (Table III). For different functionals, the prefabricated parameters are also different,86 so the dispersion energies of B3LYP-D3 is larger than that of PBE-D3 as seen in Table VI. From this section, we can see that the M06-2X is likely not suitable for Gd@C82 because its energy gap could get too large.

TABLE V. The septet-nonet gap and Gd-C distances of Gd@C82 calculated at M06-2X. ∆E = E(nonet) − E(septet). Gd-C (Å) Expt. = 2.49(3)

Method EXC M06-2X

Basis set

∆E (meV) Expt. = 1.79

Septet

Nonet

CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP

1017.63 1037.21 945.43 935.72 895.91 902.85 881.46 835.12 340.78 417.79 432.74 846.82

2.69, 2.67 2.71, 2.70 2.69, 2.68 2.71, 2.69 2.70, 2.68 2.70, 2.68 2.71, 2.70 2.70, 2.67 2.72, 2.70 2.70, 2.68 2.68, 2.67 2.71, 2.69

2.69, 2.68 2.72, 2.71 2.71, 2.69 2.72, 2.70 2.71, 2.69 2.69, 2.67 2.72, 2.70 2.71, 2.69 2.72, 2.70 2.71, 2.69 2.69, 2.67 2.71, 2.69


244306-8

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

TABLE VI. The septet-nonet gap and Gd-C distances of Gd@C82 calculated at PBE-D3 and B3LYP-D3. ∆E = E(nonet) − E(septet). Edisp is the dispersion energy in total energy. Gd’s basis set is CEP-121G.

EXC

Basis set

Septet

PBE-D3

3-21G 4-31G 4-31G* 6-31G 6-31G* cc-pvdz cc-pvtz def2-SVP 3-21G 4-31G 4-31G* 6-31G 6-31G* cc-pvdz cc-pvtz def2-SVP

−2957.31 −2959.92 −2957.57 −2958.79 −2956.87 −2956.18 −2956.05 −2957.51 −4750.07 −4755.28 −4750.67 −4753.00 −4748.69 −4746.38 −4748.95 −4749.53

B3LYP-D3

Gd-C (Å) Expt. = 2.49(3)

Edisp (meV)

Method

Nonet

∆E (meV) Expt. = 1.79

Septet

Nonet

−2957.45 −2960.09 −2957.73 −2958.96 −2957.02 −2956.35 −2956.23 −2957.68 −4750.28 −4755.58 −4750.93 −4753.30 −4748.96 −4746.66 −4749.30 −4749.83

3.10 1.97 2.17 1.44 1.71 1.34 2.38 2.07 1.02 −0.53 −0.33 −1.19 −0.90 −1.33 −0.57 −0.64

2.49, 2.49 2.48, 2.49 2.47, 2.47 2.48, 2.49 2.47, 2.47 2.46, 2.47 2.46, 2.46 2.47, 2.47 2.50, 2.50 2.49, 2.50 2.48, 2.48 2.50, 2.51 2.48, 2.49 2.48, 2.48 2.47, 2.47 2.48, 2.49

2.49, 2.49 2.48, 2.49 2.47, 2.47 2.48, 2.49 2.47, 2.47 2.46, 2.47 2.46, 2.46 2.47, 2.47 2.50, 2.50 2.49, 2.50 2.48, 2.48 2.50, 2.51 2.48, 2.48 2.48, 2.48 2.47, 2.47 2.48, 2.49

The DFT-D3 method almost does not change the conclusion on the geometric parameters and energy order compared to the corresponding standard DFT. Therefore, getting the correct ground state still depends on the appropriate standard DFT itself in the DFT-D3 method. B. The electronic structure of the Gd@C82 ground state

Based on what observed above, the calculated geometries and septet-nonet gaps of the PBE and PW91 functionals are in good agreement with the experimental results. Here, we analyzed the electronic structure of Gd@C82 at the PBE level in order to further check the reliability of pure functionals. The PW91 results are very similar to PBE as seen in Table VII. As shown in Fig. 2, the lowest unoccupied molecular orbital (LUMO) of Gd@C82 is contributed mainly by the cage (97.63%) and a small account of Gd (1.62%). Meanwhile, the highest occupied molecular orbital (HOMO) is the spin down, and the composition is similar to LUMO (cage97.99%, Gd-1.31%). Further examination on the occupied molecular orbitals (MOs) shows that no obvious 5d- or 6slike orbitals located on the Gd atom. And the hybridization between Gd and C82 are very small. As we know the electronic configuration of the isolated Gd atom is 4f 75d16s2, which indicates the 5d and 6 s electrons (3 electrons) of Gd are mostly donated to the cage suggesting the ionic interaction between Gd and C82 is important. This conclusion fits well with the Gd3+@C823− ionic model.22 The predicted HOMO-LUMO gap is about 0.10 eV. About 14 eV lower than HOMO, there are 7 single-occupied MOs located on the Gd atom which reflects the characteristic of the Gd-4f atomic orbital (the range is from –18.67 eV to –18.64 eV in Fig. 2), indicating the halffull Gd-4f shell is chemically inert and only contributes a spin magnetic moment.21 Globally, only one HOMO electron (spin down) and seven Gd-4f electrons (spin up) contribute to a septet ground state (corresponding to six net spins) for

Gd@C82. The spin density distribution is revealed in Fig. 3, which also clearly reflects the anti-ferromagnetic coupling of the net spins between Gd and the cage. These electronic structure features can well explain previous experimental observations, and the MOs information is also consistent with the Local Density Approximation plus Hubbard U (LDA+U) results.25 The orbital energy values vary from different functionals or basis sets, however the corresponding orbital energy differences from different methods are comparable. The MOs energy differences of LUMO-HOMO and HOMO-Gd(4f)1st as well as the net spins and Mulliken charges on the Gd atom are summarized in Table VII. We can see that at each C’s basis set, PBE could always get about 0.1 eV of LUMO-HOMO gap and 14 eV of HOMO-Gd(4f)1st gap. Additionally, the net spins on Gd is about 7.10 at each basis set. The Mulliken charge analysis is depended on the choice of the basis set,87,88 so the Mulliken charge is about +1.13-1.82e on Gd. These results reveal that the basis set of the carbon could have almost no effect on the electronic structure of Gd@C82 when the appropriate functional was chosen. From all the data presented thus far, it is clear that the reliability of pure functionals for Gd@C82 has been confirmed. We suggest that both PBE and PW91 functionals could be applied in relevant electronic state studies of Gd@C82 and its derivatives. In consideration of the geometry, energy gap between ground state (septet) and the next excited state (nonet), and the computational efficiency, the calculated results from medium-sized 6-31G*, etc., basis set are reasonably close to the experimental results, hence, these basis sets are also recommended. Although a previous study showed that pure-DFT overly stabilize low-spin states for iron complexes,89 we should note that both pure-GGA and hybrid-GGA predicted about +7.0 (spin up) net spins on Gd and about ±1.0 (spin up or down, corresponding to nonet or septet) net spins on the cage, hence the different spin states of Gd@C82 are not dominated by the


244306-9

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

TABLE VII. MOs energy differences and population analysis of Gd@C82 (septet). Method

MOs energy differences (eV)

Gd-populations

EXC

Basis set

LUMO-HOMO

HOMO-Gd(4f)1st a

Net spin

Charge

PBE

CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼def2-TZVP CEP121G∼3-21G CEP121G∼4-31G CEP121G∼4-31G* CEP121G∼6-31G CEP121G∼6-31G* CEP121G∼6-31+G* CEP121G∼6-311G CEP121G∼6-311G* CEP121G∼6-311+G* CEP121G∼cc-pvdz CEP121G∼cc-pvtz CEP121G∼def2-SVP CEP121G∼def2-TZVP

0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10

14.09 14.23 14.04 14.24 14.04 14.03 14.22 14.02 14.03 14.01 13.96 14.04 14.01 14.22 14.31 14.11 14.31 14.11 14.10 14.30 14.10 14.11 14.11 14.07 14.12 14.09

7.09 7.10 7.10 7.10 7.10 7.11 7.10 7.10 7.11 7.10 7.09 7.09 7.08 7.09 7.10 7.10 7.11 7.10 7.11 7.10 7.10 7.11 7.10 7.10 7.10 7.09

1.82 1.32 1.30 1.24 1.23 1.14 1.43 1.36 1.68 1.32 1.29 1.30 1.37 1.82 1.31 1.29 1.23 1.21 1.13 1.41 1.35 1.67 1.31 1.28 1.27 1.39

PW91

a The

Gd(4f)1st represents the first local 4f-like MO on Gd, such as the MO at -18.64 eV in Fig. 2.

metal Gd, but determined by the coupling between the unpaired electrons of the inner Gd and the cage. Our results show that pure-GGA could accurately reproduce the anti-ferromagnetic behavior and predict the septet ground state. Therefore, pureDFT methods present advantages in electron coupling modes predictions in Gd@C82.

FIG. 2. The molecular orbitals of the ground state (septet) of Gd@C82 calculated at PBE/CEP-121G∼6-31G*. (The 4-31G* result is very similar to 6-31G* and was presented in the supplementary material (Ref. 90).) The red and green lines represent the α-occupied and α-unoccupied MOs, respectively. The blue lines represent the β-occupied MOs. The blue and yellow areas indicate the positive and negative sign of the wave functions, respectively. Isovalue = 0.02.

C. Isomers with Gd at different positions

Previously, in both experimental and theoretical studies, the Gd atom was considered to possibly locate at different sites in the C82 cage. Therefore, we further optimized the geometries of different Gd@C82 isomers at the PBE/CEP-121G∼631G*(4-31G*) level. As seen in Fig. 4, Gd@C82-R6, in which

FIG. 3. The spin density of Gd@C82 ground state (septet). The purple and orange areas represent the upward and downward net spin electrons, respectively. (Spin density threshold value is set to 0.002.)


244306-10

Dai et al.

J. Chem. Phys. 141, 244306 (2014)

2.

3.

FIG. 4. Calculated total energies and ground states of different Gd@C82 isomers with Gd located at different sites in C82 fullerene at PBE/CEP-121G∼631G* level. (The values in parentheses are from PBE/CEP-121G∼4-31G* level.)

the Gd locates above the hexagon on the C2 axis as discussed above, is the most stable structure with a septet ground state. In contrast, the Gd@C82-B structure where the Gd locates above the C–C bond on the C2 axis is a second-order saddle point (with two imaginary frequencies) and is about 1.95 eV higher than Gd@C82-R6 in total energy. The Gd@C82-B1 can be formed by shifting the Gd away from the C2 axis relative to Gd@C82-B, which is a stationary point. However, it is about 1.28 eV higher than Gd@C82-R6. The three geometries of Gd@C82 in Fig. 4 have raised controversies in determining the most stable structure. Based on DFT calculations, we confirmed that the Gd atom tends to locate above the unique hexagon on the C2 axis of C82 to form a normal M@C82 geometry. Other geometries with Gd at different sites are either unstable or with higher total energy.

IV. CONCLUSION

In summary, calculation results of Gd@C82 are generally sensitive to the choices of calculation methods due to its small energy gap between the ground state and the next excited state. Our study showed that DFT/ECP method could reasonably describe the electronic state of Gd@C82. These methodological enlightenments can be summarized in 6 aspects: 1. The selections of functionals are decisive for the quality of results. Generally, pure functionals, such as PBE and PW91, could determine the correct septet ground state. Both the geometric parameters and septet-nonet gap from PBE and PW91 agree with the experimental results. Therefore, these pure functionals are recommended in the relevant electronic structures studies of Gd@C82. Here, it should also be noticed that pure-DFT methods suffer from the unavoidable self-interaction error. But more importantly, the accurate prediction of the ground states is the basic precondition for studying other properties of EMFs. Since pure-DFT have better performance than hybrid-DFT in this aspect as supported by our results on Gd@C82, pure-

4.

5.

6.

DFT methods are preferentially recommended in electronic state predictions of Gd-based EMFs. When the pure functionals are used, the common mediumsized double-split 6-31G*, etc., are good enough to calculate reasonable septet-nonet gap and Gd-C distance. Enlarging the basis set or introducing the diffuse functions do not improve the results (actually making them slightly worse). The introduction of the polarization functions will decrease the Gd-C distance by 0.01–0.03 Å and increase the energy gap by 0.2–0.3 meV, whereas the diffuse functions will increase the Gd-C distance by about 0.01 Å and decrease the energy gap by about 0.5 meV. The ground states from hybrid functionals strongly depend on C’s basis sets. Only when small basis sets are selected, the hybrid functionals could obtain qualitatively correct septet ground state. Enlarging C’s basis sets may not necessarily lead to reasonable results. Similar phenomenon was also found in our previous study.49 LSDA could determine the septet ground state with less dependence on C’s basis sets, however LSDA overestimates the septet-nonet energy gap and underestimates the Gd-C distance. M06-2X predicts poor septet-nonet energy gap and overestimates the Gd-C distance. For double-split basis sets, the ZPE correction for septet is larger than that for nonet and results in 0.4–0.6 meV reduction on the septet-nonet energy gap. For triplet-split 6311G basis sets, the ZPE corrections reduce the energy gap by more than 1 meV, and for the 6-311G*, the reduction is about 0.3–0.4 meV. Considering that the total septet-nonet energy gap is only 1.79 meV according to the experiment, the ZPE contribution cannot be neglected. For the DFT-D3 method, the introduction of dispersion effect almost does not change the geometric parameters as compared to the standard DFT. The dispersion energy brings a 0.15 meV reduction to the septet-nonet energy gap due to the slight difference in the geometries between the septet and nonet states.

In short, determining reliable calculation methods is an important issue in f-element EMFs area. The EMF Gd@C82 turns out to be an ideal case for benchmarking tests in both theoretical and experimental studies due to its very small energy gap, which has contributed to the difficulty of determining the ground electronic state and most stable geometry. Our conclusions could provide important references for subsequent studies of Gd@C82 and its derivatives, as well as the structural and electronic features of new f-element EMFs.

ACKNOWLEDGMENTS

We would like to thank the support of the National Science Foundation of China (11374004, 21320102003) and Science and Technology Development Program of Jilin Province of China (20150519021JH). Z.W. acknowledges the assistance of the Fok Ying Tung Education Foundation (142001) and the High Performance Computing Center (HPCC) of Jilin University.


244306-11 1A.

Dai et al.

A. Popov, S. Yang, and L. Dunsch, Chem. Rev. 113(8), 5989–6113 (2013). 2H. Shinohara, Rep. Prog. Phys. 63(6), 843–892 (2000). 3A. G. Starikov, O. A. Gapurenko, A. L. Buchachenko, A. A. Levin, and N. N. Breslavskaya, Russ. J. Gen. Chem. 78(4), 793–810 (2008). 4A. L. Kutepov, J. Alloys Compd. 444-445, 174–176 (2007). 5K. I. M. Ingram, M. J. Tassell, A. J. Gaunt, and N. Kaltsoyannis, Inorg. Chem. 47(17), 7824–7833 (2008). 6P. Wåhlin, C. Danilo, V. Vallet, F. Réal, J.-P. Flament, and U. Wahlgren, J. Chem. Theory Comput. 4(4), 569–577 (2008). 7P. Pyykko, Chem. Rev. 88(3), 563–594 (1988). 8M. Dolg, in Theoretical and Computational Chemistry, edited by S. Peter (Elsevier, 2002), Vol. 11, pp. 793–862. 9L. Gagliardi, P. Pyykko, and B. O. Roos, Phys. Chem. Chem. Phys. 7(12), 2415–2417 (2005). 10L. Gagliardi and B. O. Roos, Nature 433(7028), 848–851 (2005). 11M. F. Zalazar, V. M. Rayon, and A. Largo, J. Phys. Chem. A 116(11), 2972–2977 (2012). 12Y. Gao, X. Dai, S.-G. Kang, C. A. Jimenez-Cruz, M. Xin, Y. Meng, J. Han, Z. Wang, and R. Zhou, Sci. Rep. 4, 5862 (2014). 13B. O. Roos and P.-A. Malmqvist, Phys. Chem. Chem. Phys. 6(11), 2919 (2004). 14X. Wu and X. Lu, J. Am. Chem. Soc. 129(7), 2171–2177 (2007). 15C. Chen, G. Xing, J. Wang, Y. Zhao, B. Li, J. Tang, G. Jia, T. Wang, J. Sun, L. Xing, H. Yuan, Y. Gao, H. Meng, Z. Chen, F. Zhao, Z. Chai, and X. Fang, Nano Lett. 5(10), 2050–2057 (2005). 16S.-G. Kang, T. Huynh, and R. Zhou, Sci. Rep. 2, 957 (2012). 17S.-G. Kang, G. Zhou, P. Yang, Y. Liu, B. Sun, T. Huynh, H. Meng, L. Zhao, G. Xing, C. Chen, Y. Zhao, and R. Zhou, Proc. Natl. Acad. Sci. U. S. A. 109(38), 15431–15436 (2012). 18J. Meng, X. Liang, X. Chen, and Y. Zhao, Integr. Biol. 5(1), 43–47 (2013). 19H. Funasaka, K. Sakurai, Y. Oda, K. Yamamoto, and T. Takahashi, Chem. Phys. Lett. 232(3), 273–277 (1995). 20H. Funasaka, K. Sugiyama, K. Yamamoto, and T. Takahashi, J. Phys. Chem. 99(7), 1826–1830 (1995). 21T. Suzuki, K. Kikuchi, F. Oguri, Y. Nakao, S. Suzuki, Y. Achiba, K. Yamamoto, H. Funasaka, and T. Takahashi, Tetrahedron 52(14), 4973–4982 (1996). 22K. Suenaga, S. Iijima, H. Kato, and H. Shinohara, Phys. Rev. B 62(3), 1627–1630 (2000). 23H. J. Huang, S. H. Yang, and X. X. Zhang, J. Phys. Chem. B 103(29), 5928–5932 (1999). 24H. Huang, S. Yang, and X. Zhang, J. Phys. Chem. B 104(7), 1473–1482 (2000). 25A. Sebetci and M. Richter, J. Phys. Chem. C 114(1), 15–19 (2009). 26K. Furukawa, S. Okubo, H. Kato, H. Shinohara, and T. Kato, J. Phys. Chem. A 107(50), 10933–10937 (2003). 27H. Giefers, F. Nessel, S. I. Györy, M. Strecker, G. Wortmann, Y. S. Grushko, E. G. Alekseev, and V. S. Kozlov, Carbon 37(5), 721–725 (1999). 28E. Nishibori, M. Takata, M. Sakata, H. Tanaka, M. Hasegawa, and H. Shinohara, Chem. Phys. Lett. 330(5–6), 497–502 (2000). 29E. Nishibori, M. Takata, M. Sakata, M. Inakuma, and H. Shinohara, Chem. Phys. Lett. 298(1–3), 79–84 (1998). 30B. Sun, M. Li, H. Luo, Z. Shi, and Z. Gu, Electrochim. Acta 47(21), 3545–3549 (2002). 31E. Nishibori, K. Iwata, M. Sakata, M. Takata, H. Tanaka, H. Kato, and H. Shinohara, Phys. Rev. B 69(11), 113412 (2004). 32L. Liu, B. Gao, W. Chu, D. Chen, T. Hu, C. Wang, L. Dunsch, A. Marcelli, Y. Luo, and Z. Wu, Chem. Commun. 2008(4), 474–476. 33T. Akasaka, T. Kono, Y. Takematsu, H. Nikawa, T. Nakahodo, T. Wakahara, M. O. Ishitsuka, T. Tsuchiya, Y. Maeda, M. T. H. Liu, K. Yoza, T. Kato, K. Yamamoto, N. Mizorogi, Z. Slanina, and S. Nagase, J. Am. Chem. Soc. 130(39), 12840–12841 (2008). 34M. Suzuki, X. Lu, S. Sato, H. Nikawa, N. Mizorogi, Z. Slanina, T. Tsuchiya, S. Nagase, and T. Akasaka, Inorg. Chem. 51(9), 5270–5273 (2012). 35K. Kobayashi and S. Nagase, Chem. Phys. Lett. 282(3–4), 325–329 (1998). 36L. Senapati, J. Schrier, and K. B. Whaley, Nano Lett. 4(11), 2073–2078 (2004). 37L. Senapati, J. Schrier, and K. B. Whaley, Nano Lett. 5(11), 2341–2341 (2005). 38L. Wang and D. Yang, Nano Lett. 5(11), 2340–2340 (2005). 39N. Mizorogi and S. Nagase, Chem. Phys. Lett. 431(1–3), 110–112 (2006). 40A. A. Popov and L. Dunsch, Chem. Eur. J. 15(38), 9707–9729 (2009).

J. Chem. Phys. 141, 244306 (2014) 41M.

J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. Montgomery, Jr., J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, T. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, O. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski, and D. J. Fox, Gaussian’09, Revision D.01 (Gaussian, Inc., Wallingford CT, 2013). 42T. R. Cundari and W. J. Stevens, J. Chem. Phys. 98(7), 5555–5565 (1993). 43J. G. Rodríguez-Zavala, I. Padilla-Osuna, and R. A. Guirado-López, Phys. Rev. B 78(15), 155426 (2008). 44M. Dolg, H. Stoll, and H. Preuss, J. Chem. Phys. 90(3), 1730–1734 (1989). 45X. Cao and M. Dolg, J. Chem. Phys. 115(16), 7348–7355 (2001). 46J.-P. Dognon, C. Clavaguéra, and P. Pyykkö, J. Am. Chem. Soc. 131(1), 238–243 (2008). 47X. Dai, C. Cheng, W. Zhang, M. Xin, P. Huai, R. Zhang, and Z. Wang, Sci. Rep. 3, 1341 (2013). 48J. Han, X. Dai, C. Cheng, M. Xin, Z. Wang, P. Huai, and R. Zhang, J. Phys. Chem. C 117(50), 26849–26857 (2013). 49M. Xin, X. Dai, B. Huang, Y. Meng, W. Feng, M. Jin, Z. Wang, and R.-Q. Zhang, Chem. Phys. Lett. 585, 107–111 (2013). 50M. Xin, X. Dai, J. Han, M. Jin, C. A. Jimenez-Cruz, D. Ding, Z. Wang, and R. Zhou, RSC Adv. 4(57), 30074–30080 (2014). 51X. Dai, M. Xin, Y. Meng, J. Han, Y. Gao, W. Zhang, M. Jin, Z. Wang, and R.-Q. Zhang, Carbon 78, 19–25 (2014). 52X. Dai, J. Han, Y. Gao, and Z. Wang, ChemPhysChem 15, 3871–3876 (2014). 53J. S. Binkley, J. A. Pople, and W. J. Hehre, J. Am. Chem. Soc. 102(3), 939–947 (1980). 54R. Ditchfield, W. J. Hehre, and J. A. Pople, J. Chem. Phys. 54(2), 724–728 (1971). 55W. J. Hehre, R. Ditchfield, and J. A. Pople, J. Chem. Phys. 56(5), 2257–2261 (1972). 56P. C. Hariharan and J. A. Pople, Theor. Chem. Acc. 28, 213–222 (1973). 57T. Clark, J. Chandrasekhar, G. W. Spitznagel, and P. V. R. Schleyer, J. Comp. Chem. 4, 294 (1983). 58R. Krishnan, J. S. Binkley, R. Seeger, and J. A. Pople, J. Chem. Phys. 72(1), 650–654 (1980). 59F. Weigend, Phys. Chem. Chem. Phys. 8(9), 1057–1065 (2006). 60F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys. 7(18), 3297–3305 (2005). 61T. H. Dunning, J. Chem. Phys. 90(2), 1007–1023 (1989). 62R. A. Kendall, T. H. Dunning, and R. J. Harrison, J. Chem. Phys. 96(9), 6796–6806 (1992). 63J. C. Slater, Quantum Theory of Molecules and Solids: The Self-Consistent Field for Molecules and Solids (McGraw-Hill, 1963). 64P. Hohenberg and W. Kohn, Phys. Rev. 136(3B), B864–B871 (1964). 65W. Kohn and L. J. Sham, Phys. Rev. 140(4A), A1133–A1138 (1965). 66S. H. Vosko, L. Wilk, and M. Nusair, Can. J. Phys. 58(8), 1200–1211 (1980). 67J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77(18), 3865–3868 (1996). 68J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 78(7), 1396–1396 (1997). 69J. P. Perdew, in Electronic Structure of Solids, edited by P. Ziesche and H. Eschrig (Akademie Verlag, Berlin, 1991), Vol. 11. 70J. P. Perdew, K. Burke, and Y. Wang, Phys. Rev. B 54(23), 16533–16539 (1996). 71J. P. Perdew, J. A. Chevary, S. H. Vosko, K. A. Jackson, M. R. Pederson, D. J. Singh, and C. Fiolhais, Phys. Rev. B 46(11), 6671–6687 (1992). 72J. P. Perdew, J. A. Chevary, S. H. Vosko, K. A. Jackson, M. R. Pederson, D. J. Singh, and C. Fiolhais, Phys. Rev. B 48(7), 4978–4978 (1993). 73K. Burke, J. P. Perdew, and Y. Wang, in Electronic Density Functional Theory: Recent Progress and New Directions, edited by J. F. Dobson, G. Vignale, and M. P. Das (Plenum, 1998). 74A. D. Becke, Phys. Rev. A 38(6), 3098–3100 (1988). 75C. Lee, W. Yang, and R. G. Parr, Phys. Rev. B 37(2), 785–789 (1988). 76B. Miehlich, A. Savin, H. Stoll, and H. Preuss, Chem. Phys. Lett. 157(3), 200–206 (1989).


244306-12 77J.

Dai et al.

P. Perdew, Phys. Rev. B 33(12), 8822–8824 (1986). D. Becke, Chem. Phys. 98(7), 5648–5652 (1993). 79C. Adamo and V. Barone, J. Chem. Phys. 110(13), 6158–6170 (1999). 80S. Grimme, J. Antony, S. Ehrlich, and H. Krieg, J. Chem. Phys. 132(15), 154104 (2010). 81Y. Zhao and D. G. Truhlar, Theor. Chem. Acc. 120(1–3), 215–241 (2007). 82A. R. Khamatgalimov and V. I. Kovalenko, J. Phys. Chem. A 115(44), 12315–12320 (2011). 83S. I. Gorelsky, J. Chem. Theory Comput. 8(3), 908–914 (2012). 84X. Liu, L. Li, B. Liu, D. Wang, Y. Zhao, and X. Gao, J. Phys. Chem. A 116(47), 11651–11655 (2012). 78A.

J. Chem. Phys. 141, 244306 (2014) 85A.

A. Popov, S. M. Avdoshenko, A. M. Pendas, and L. Dunsch, Chem. Commun. 48(65), 8031–8050 (2012). 86See http://www.gaussian.com/g_tech/g_ur/k_dft.htm for Gaussian 09 user’s reference. 87P. Belanzoni, E. van Lenthe, and E. J. Baerends, J. Chem. Phys. 114(10), 4421–4433 (2001). 88A. E. Reed, R. B. Weinstock, and F. Weinhold, J. Chem. Phys. 83(2), 735–746 (1985). 89M. Swart, A. R. Groenhof, A. W. Ehlers, and K. Lammertsma, J. Phys. Chem. A 108(25), 5479–5483 (2004). 90See supplementary material at http://dx.doi.org/10.1063/1.4904389 for DFT/ANO/ECP, DKH-DFT, and ZORA-DFT calculations.


The electronic ground state of [Fe(CO)3 (NO)](-) : a spectroscopic and theoretical study.

Ground-state densities from the Rayleigh-Ritz variation principle and from density-functional theory.

Simple Fully Nonlocal Density Functionals for Electronic Repulsion Energy.

Communication: The ground electronic state of Si2C: Rovibrational level structure, quantum monodromy, and astrophysical implications.

Theoretical study of the first electronic state of the FCO2 radical.

Theoretical investigation of the effects of doping on the electronic structure and thermoelectric properties of ZnO nanowires.

Guanidine and guanidinium cation in the excited state--theoretical investigation.

3-Center-5-electron boryl radicals with σ(0)π(1) ground state electronic structure.

Thermodynamic ground state of MgB6 predicted from first principles structure search methods.

Investigation of VO-salophen complexes electronic structure.

Experimental and theoretical investigation on the outer valence electronic structure of cyclopropylamine by (e, 2e) electron momentum spectroscopy.

Theoretical and Experimental Photoelectron Spectroscopy Characterization of the Ground State of Thymine Cation.

First-principles study of ground-state properties of U2Mo.

Theoretical and Experimental Studies on the Crystal Structure, Electronic Structure and Optical Properties of SmTaO₄.

Systematic theoretical investigation of the zero-field splitting in Gd(III) complexes: wave function and density functional approaches.

Density-functional theory of the electronic structure of molecules.

Synthesis and electronic structure of the first cyaphide-alkynyl complexes.

Transition-state searches in metal clusters by first-principle methods.

Experimental and theoretical investigation of the spectroscopic and electronic properties of pyrazolyl ligands.

Deviations from piecewise linearity in the solid-state limit with approximate density functionals.

Electronic structure investigation of the evanescent AtO(+) ion.

First principle investigation of the linker length effects on the thermodynamics of divalent pseudorotaxanes.

A theoretical investigation of mixing thermodynamics, age-hardening potential, and electronic structure of ternary M(1)1-x M(2)xB2 alloys with AlB2 type structure.

Density-based partitioning methods for ground-state molecular calculations.