Article pubs.acs.org/JPCC
Performance of Dispersion-Corrected DFT for the Weak Interaction between Aromatic Molecules and Extended Carbon-Based Systems Doreen Mollenhauer,*,†,‡ Claudia Brieger,‡ Elena Voloshina,§ and Beate Paulus‡ †
Institute of Physical Chemistry, Justus-Liebig-University Giessen, Heinrich-Buff-Ring 58, 35392 Giessen, Germany Institute for Chemistry and Biochemistry, Freie Universität Berlin, Takustr. 3, 14195 Berlin, Germany § Institut for Chemistry, Humboldt-Universität zu Berlin, Unter den Linden 6, 10099 Berlin, Germany ‡
S Supporting Information *
ABSTRACT: The performance of different DFT approaches in combination with dispersion correction is studied for the interaction between aromatic molecules and extended carbonbased materials on the example of the pyridine−graphene system. The basic interaction is modeled using graphene fragments of increasing size as well as periodic boundary conditions. Different DFT-D2/D3 methods are tested for small and medium fragment systems in comparison to wavefunction-based CCSD(T) and SCS-MP2 approaches. Furthermore, the adsorption energy between pyridine and extended graphene sheets or periodic modeled graphene calculated by DFT-D2/D3 or nonlocal correlation functionals (vdW-DF) is compared to experimental values. The study of DFT-D performance along different scales reveals the dispersion correction as too strong along increasing graphene fragment sizes. Finally, this leads to different methodology advice for small and extended pyridine−graphene systems.
■
INTRODUCTION Graphene has been found to show unique electronic, mechanical, and transport properties which are of high potential for applications in the fields of electrochemical sensing, nanoelectronic device materials, environmental filters, drug delivery, or energy technologies.1−4 Considering aromatic molecule binding modes to graphene, the interaction of the building block pyridine to graphene is of special interest from the medical applications point of view, such as diagnostic devices or drug-delivery systems. Here, graphene can be used for the passivation of metal nanoparticles in order to improve stability or avoid toxicity effects (e.g., ref 5). Further, the pyridine−graphene interaction has importance in pyridinecapping CdSe−graphene pholtovoltaic devices, nitrogen-doped graphene for sensors, or Li-ion batteries.6−8 In our recent theoretical study using the PBE-D2 method we found the adsorption of pyridine to graphene and graphite is preferable in a flat orientation9 in contrast to the normal mode interaction of pyridine to several noble metallic surfaces.10−14 Furthermore, the adsorption energy for the pyridine−graphene system was significantly weaker than for pyridine to the metal surfaces.10−14 The main contribution to the adsorption energy (about 90%) comes from the London dispersion (C6R−6); the band structure of graphene was determined to be unchanged by this weak interaction. Other theoretical investigations regarding pyridine-like molecules in interaction with graphene are aminotriazines on graphene using a PBE-D2 approach.15 The authors found an enhancement of the adsorption energy due to © 2015 American Chemical Society
NR2 groups on benzene, pyridine, or triazine cores. Pyridinelike supramolecular building blocks were studied also on graphite. Here the DFT-D2 and the DFT-D3 method agreed reasonably well with experimental values, whereas force-field approaches varied strongly in determined adsorption energies.16 Another investigation treats the pyridine adsorption to zigzag graphene nanoribbons; here also the PBE-D2 approach leads to increasing adsorption energies and shorter pyridine graphene distances compared to pure PBE.17 Furthermore, nicotin was studied on graphene using a hybrid functional, as well as nucleobase graphene systems using different DFT-D2/ D3, vdW-DF, and TS approaches and SCS-MP2.18−20 In all cases the dispersion correction plays a crucial role for the appropriate description of the adsorption energy. Due to the high importance of the dispersion for pyridinelike molecules to graphene the aim of this paper is to investigate the accuracy of dispersion correction methods over different scales and to compare the outcome on the one hand to benchmarking methods and on the other hand to experiments. In this work we investigate the bonding of a single pyridine molecule to the graphene surface. In a first step the graphene surface is modeled as finite fragment with increasing size. Dangling bonds are saturated with H atoms. We use the highly symmetric parallel on-top arrangement (Figure Received: November 12, 2014 Revised: January 2, 2015 Published: January 20, 2015 1898
DOI: 10.1021/jp5113312 J. Phys. Chem. C 2015, 119, 1898−1904
Article
The Journal of Physical Chemistry C
Figure 1. Used pyridine on-top graphene fragment structures, sized C6H6 (a), C24H12 (b), C54H18(c), C96H25 (d), and C105H30 (e).
1) as preferred flat orientation for the dispersion-dominated interaction. This adsorption mode can be integrated much better in our aim of investigation; however, the hollow adsorption site was found to be energetically more preferred (by 16%) in previous investigations.9 The bonding of pyridine to the medium-size graphene fragment C24H12 serves to test different DFT approaches (LDA, GGA, hybrid, and double hybrid functionals) in conjunction with the semiempirical dispersion correction scheme of Grimme (D2 and D3). The results are compared to wave function based methods MP2 and SCS-MP2 to determine the accuracy of the functionals and dispersion correction. To prove the accuracy of the (SCS)-MP2 calculations, the highly correlated wave-function-based CCSD(T) method has been applied to the C 6H 6 −pyridine interaction. Thereby also the basis set was varied, and an extrapolation was performed to the complete basis set limit. Second, the size of the saturated graphene fragment is varied with C6H6, C24H12, C54H18, C96H24, and C150H30 (Figure 1, a− e) to show the convergence of the adsorption energy for different functional/dispersion correction combinations and compare to periodic results with different vdW-DF, vdW-DF2, and “opt” functionals. The approach allows us to determine statements by methodological and experimental comparison of the scaling effects of the dispersion correction scheme for aromatic molecules in interaction to graphene.
functional B2PLYP42 are employed. Medium- and long-range dispersion effects are added within D2 and D3 schemes to the total energies for all functionals and noted as −D2 or −D3 at the functional when applied.22 Furthermore, we varied the damping function (zero, BJ) which determines the short- and medium-range dispersion correction in order to avoid near singularities for small distances.43 On one hand, the damping to zero for short ranges is utilized, and on the other hand the damping proposed by Becke and Johnson (BJ) is applied. As the inclusion of three-body nonadditive dispersion terms can be particularly important for larger systems we have tested this and marked it by E3 if applied. In the calculations the CN and CC distances of pyridine as well as the graphene fragment are kept fixed at experimental values.44 Previous calculations indicated just very small structural changes (less than 1%) for pyridine as well as graphene due to adsorption.9 The pyridine−graphene distance (difference between pyridine and graphene plane) and adsorption energy between pyridine and graphene fragments are calculated by potential energy curves, which are fitted to polynomials of sixth order around the minimum region. The basis set superposition error (BSSE) has been corrected by a counterpoise correction for all pyridine graphene fragment calculations.45 For MP2 and SCS-MP2 basis set extrapolation of the correlation energy is done by the formula of ref 45. The periodic DFT calculations are carried out with the VASP program46 using the projector augmented wave method47 and plane wave basis set. The plane wave kinetic energy cutoff is set to 500 eV. The calculations are performed using a (4 × 4) surface periodicity that corresponds to an adsorbate concentration of approximately 1.2 × 10−2 Å−2. In the total energy calculations the k-meshes for sampling of the supercell Brillouin zone are chosen to be as dense as 36 × 36 k-mesh, when folded up to the simple graphene unit cell. The vacuum region between graphene layers is set to ca. 24 Å. The relevance of long-range van der Waals interactions on the molecule−surface adsorption process is investigated by means of (i) a semiempirical DFT-D2 approach proposed by Grimme or (ii) nonlocal correlation functionals: original vdW-DF,48 its three modifications, where the exchange functional is optimized for the dispersion part (optPBE-vdW, optB88-vdW, and optB86b-vdW),49 and the vdW-DF2 functional.50
■
METHOLOGY All nonperiodic calculations are performed with the program package MOLPRO,21 and the dispersion correction is applied by using the dftd3 code22 from the Grimme group. The Dunning-type correlation consistent valence basis sets (ccpVXZ) with X = D, T, Q are used for the elements H, C, and N.23,24 As a wave function-based approach, the Hartree−Fock method and the second-order Møller−Plesset perturbation method (MP2) are utilized.25,26 Additionally, the spincomponent-scaled second-order Møller−Plesset perturbation method (SCS-MP2) is applied.27 It was shown that this method gives good agreement to CCSD(T).28−32 The MP2 approach can be applied for the system as fragments described with a finite basis set are used. For the larger pyridine−graphene systems, density fitting is applied for the MP2 approaches.33 Test calculations reveal differences in adsorption energy smaller than 1 meV due to this approximation. Also, different approximations of DFT, namely, local density approximation PWLDA,34 generalized gradient approximation (GGA), here PBE,35 BP86,36,37 and B9738 functionals, metaGGA TPSS,39 hybrid functional B3LYP,40,41 and double hybrid
■
RESULTS AND DISCUSSION Wave-Function-Based Approaches for the Finite Pyridine−Graphene Fragment System. Wave-functionbased methods have served to describe the pyridine−C6H6, 1899
DOI: 10.1021/jp5113312 J. Phys. Chem. C 2015, 119, 1898−1904
Article
The Journal of Physical Chemistry C
Table 1. Adsorption Energies (EAds in meV) and Pyridine−Graphene Distance (d0 in Å) for the Pyridine−C6H6, −C24H12, and −C54H18 Systems Calculated with the MP2 and SCS-MP2 Method and Different Basis Sets and Extrapolation to the Complete Basis Set Limit cc-pVDZ method C6H6 CCSD(T) MP2 SCS-MP2 C24H12 MP2 SCS-MP2 C54H18 MP2 SCS-MP2
cc-pVTZ
EP (DZ/TZ)
cc-pVQZ
EP (TZ/QZ)
EAds
d0
EAds
d0
EAds
d0
EAds
d0
EAds
d0
22 64 31
4.20 3.85 4.15
60 128 72
3.93 3.70 3.88
76 156 89
3.87 3.65 3.85
153 87
3.65 3.83
170 98
3.61 3.78
338 209
3.34 3.55
522 331
3.25 3.45
605 387
3.25 3.35
-
-
-
-
442 282
3.345 3.467
-
-
-
-
-
-
-
-
Figure 2. Basis set dependency and extrapolation for the pyridine−C6H6 system (a) using the CCSD(T) (b) and SCS-MP2 (c) method, respectively.
−C24H12, and −C54H18 fragments, and the outcome is given in Table 1 and Figure 2. The smallest model system, namely, C6H6−pyridine (Figure 2a), has served for CCSD(T) calculations and the proof of accuracy for the (SCS)-MP2 method [(SCS)-MP2 is used if both approaches MP2 and SCSMP2 are meant]. The interaction is very weak, but nevertheless the SCS-MP2 result is in very good agreement with the CCSD(T) one (differences smaller than 10 meV). In contrast, MP2 overestimates the interaction by a factor of 2−3. Thereby the difference in MP2 and SCS-MP2 results decreases slightly by increasing the basis set (from a factor of 2.3 to 1.7). Also, MP2 tends to overestimate the pyridine−graphene distance by about 0.2 Å compared to SCS-MP2. In all cases the HF calculations result in a nonbonding situation. Electron correlation effects determine the dispersiondominated binding situation. The increase of the basis set is very important for CCSD(T) and (SCS-)MP2 calculations, and there is a doubling of adsorption energy from applying DZ to TZ quality. However, the estimated extrapolation DZ/TZ agrees very well with the TZ/QZ one, which confirms the DZ/ TZ extrapolated values close to the basis set limit. On the basis of this result, we estimate the (SCS)-MP2 with DZ/TZ extrapolated values for the C24H12−pyridine close to the basis set limit, too. The difference between SCS-MP2 and MP2 for the C24H12−pyridine system decreases to a factor of 1.6, which remains for the C54H18−pyridine system; here the influence of the basis set and system size is very small and within the factor
of 1.6. The basis set dependency is visualized in Figure 2b,c. Augmented functions are also tested for the C6H6−pyridine system and (SCS-)MP2 with an AVDZ/AVTZ extrapolation and yield just 10 meV difference to the DZ/TZ extrapolated value. By considering the size extension of the graphene fragment from C6H6 to C24H12, the adsorption energy increases strongly by a factor of 5 and 7 for MP2 and SCS-MP2 (cc-pVDZ basis set), respectively. However, there is some contribution of basis set dependency because for triple-ζ quality the factors increase just by 4 and 5. The change in adsorption energy is much smaller for considering the extension from C24H12 to C54H18 (factor 1.3). In the following we discuss a situation close to convergence in system size regarding the adsorption energy for the pyridine−C54H18 system. In summary, the C6H6−pyridine interaction is investigated at a highly correlated wave function based level of theory, namely, CCSD(T). The results are in very good agreement with SCSMP2 results, thus this method serves for a benchmark method of larger systems. The SCS-MP2 method was already identified as an accurate method for the benzene−benzene interaction and is much improved compared to MP2.32 The DZ/TZ extrapolation provides reasonable results close to the basis set limit. The adsorption energy for the pyridine−graphene fragment interaction increases strongly from C6H6 to C24H12 and much less from C24H12 to C54H18. 1900
DOI: 10.1021/jp5113312 J. Phys. Chem. C 2015, 119, 1898−1904
Article
The Journal of Physical Chemistry C Method Comparison for the Pyridine−C24H12 System. Pyridine on C24H12 is used as a reasonable large system to determine the DFT-D2/D3 performance compared to SCSMP2 results. The dispersion correction term represents the dominant contribution to the adsorption energy when using the DFT-(D2/D3) method. All results are summarized in Figure 3
D3(BJ) regarding the adsorption energy and for B3LYP-D2, BP86-D3(zero,E3), and B2PLYP-D2 regarding the bond distance. Overall BP86-D3(zero,E3), B97-D3(BJ,E3), and B3LYP-D3(BJ) give the best agreement to SCS-MP2 considering both values in combination. Variation of the Graphene Fragment Size up to Periodic Boundary Conditions. Pyridine is investigated on graphene fragments of different size, namely, C6H6, C24H12, C54H18, C96H24, and C150H30, to reach convergence of the adsorption energy. Furthermore, periodic calculations of the single pyridine adsorption to a graphene layer have been performed. Therefore, the best working GGA functional (for pyridine−C24H12), BP86, and hybrid functional, B3LYP, are chosen as well as the GGA functional PBE because this is the most used one regarding pyridine/graphene and pyridine-like doped graphene theoretical investigations.15−17 The double-ζ basis set is used for the different calculations, and the difference to triple-ζ quality is small for the pyridine−C24H12 system within the range from 15 to 22 meV (depending on the density functional). The results are given in the Supporting Information Tables S2 and S3 and are visualized in Figure S2 and Figure 4.
Figure 3. DFT-D results and SCS-MP2 benchmark of the adsorption energy and pyridine−graphene distance regarding the pyridine− C24H12 system. The different dispersion corrections are marked and grouped together to analyze the overall trend.
and in the Supporting Information Table S1. The binding energies evaluated with the different functionals without the D2/D3 correction yield a nonbinding situation, and just PWLDA and B2PLYP result in a binding with strongly underestimated values compared to SCS-MP2. The addition of the dispersion correction to the DFT results improves the outcome significantly. Thereby, all dispersion-corrected functionals with zero-damping can be clustered to a specific range of adsorption energy and pyridine−graphene distance; further, there is a small shift of the energy values by adding the threebody nonadditive dispersion term. In contrast, the results for the BJ-damping vary much more in final binding energies as well as pyridine−graphene distance. The D2-damping shows an intermediate situation. The results are plotted in Figure 3 and marked regarding the dispersion correction. However, for the discussed trends the dispersion-corrected BP86 is an exception, and the outcome cannot be grouped to the other values. The adsorption energy is for most functionals in combination with dispersion correction best descripted using the BJdamping, therefore, without the three-body nonadditive dispersion term for B2PLYP, B3LYP, TPSS, and PBE and with this term for B97-D. The functionals B2PLYP, B3LYP, and B97-D yield variations of the adsorption energy within 40 meV (10% of the adsorption energy) regarding the SCS-MP2 benchmark. Also the BP86 functional results in very good agreement to SCS-MP2, but with zero-damping. The D3correction gives in all cases improved results in comparison to the older D2-correction. However, in contrast the D2correction results in better agreement of pyridine−graphene distance compared to D3 (an exception represents again BP86). Considering the shape of the potential energy curve using the dispersion correction with the best agreement to SCS-MP2, all curves match the benchmark curve well except the dispersion-corrected BP86 one (see in the Supporting Information Figure S1). In summary, the best agreement to benchmark results is found for B97-D3(BJ,E3), BP86-D3(zero,E3), and B3LYP-
Figure 4. Different periodic approaches (for BP86 the fragment C150H30 approach is used) for the calculation of the pyridine− graphene interaction. The gray line represents the experimental value and the gray dotted line the adjusted experimental value (for details see text).
The size-extended DFT calculations show convergence (
