Structure and binding in crystals of cage-like molecules:
hexamine and platonic hydrocarbons
Kristian Berland and Per Hyldgaard
arXiv:1010.1487v1 [physics.chem-ph] 7 Oct 2010
Department of Microtechnology and Nanoscience, MC2, Chalmers University of Technology, SE-41296 Göteborg, Sweden
(Dated: October 8, 2010)
In this paper, we show that first-principle calculations using a van der Waals density functional
(vdW-DF), [Phys. Rev. Lett. 92, 246401 (2004)] permits determination of molecular crystal
structure. We study the crystal structures of hexamine and the platonic hydrocarbons (cubane and
dodecahedrane). The calculated lattice parameters and cohesion energy agree well with experiments.
Further, we examine the asymptotic accounts of the van der Waals forces by comparing full vdW-DF
with asymptotic atom-based pair potentials extracted from vdW-DF. The character of the binding
differ in the two cases, with vdW-DF giving a significant enhancement at intermediate and relevant
binding separations. We analyze consequences of this result for methods such as DFT-D, and
question DFT-D’s transferability over the full range of separations.
PACS numbers: 61.50.Lt, 61.66.Hq, 71.15.Mb
I.
INTRODUCTION
Understanding supramolecular structure and interactions is essential for understanding many biological
processes.1 Biological and other supramolecular complexes as polymers and overlayers are sparse matter, that
is, they contain low electronic density in essential regions.
The general lack of order in these systems prohibits precise measurements of atomic structure and therefore challenges development of theoretical methods. In turn, this
makes transferable first principles schemes attractive; an
accurate account of simple periodic structures permits
accurate characterization to be made and reliable conclusions to be drawn for more complex (non-periodic)
systems. Molecular crystals and simple polymer crystals2
are ideal testing grounds for applications of first-principle
descriptions of sparse matter. Unlike most sparse matter, they constitute ordered systems and therefore leads
to unambiguous comparison between theory and experiments.
The dispersion forces underpin the cohesion of sparse
matter. Modeling of sparse matter at the electronic level,
therefore, requires that we take these effects into account.
Traditional implementations of Density-Functional theory (DFT), which parameterizes the density functional
within the local density approximation (LDA)3 or the
generalized gradient approximation (GGA)4,5 have allowed for routine modeling of matter with dense electron
distributions. However, the lack of non-local correlation
and hence dispersion forces in these implementations,
has effectively inhibited widespread use for general sparse
matter systems. Computational methods able to predict
molecular crystal structure and stability are also of value
in the pharmaceutical industry, for example in mapping
out competing crystalline phases for the fabrication of
pharmaceuticals.6 The electronic structure and response
are often required to understand and compute properties
of technological significance. Consequentially, a range of
developments has aimed to extend density functional approximations with an account of dispersive interactions
and thus capability to address sparse matter challenges.
Such approaches include non-local functionals7–12 , intermolecular perturbation theory13,14 , and the addition of
a semi-emperical atom-based account of the dispersion
forces (DFT-D).15–18 DFT-D has also been used to characterize molecular crystals.19,20
We study molecular crystal using a van der Waals density functional (vdW-DF).10 The correlation part of this
functional depends non-locally on the density and accounts for dispersion forces. vdW-DF has been used
in a series of first-principle studies of sparse matter,21
yielding for example new insight to the twist of DNA,22
nanotube bundles,23 water hexamers,24 metal-organic
frameworks,25 and the binding mechanisms at organicmetal interfaces.26–29
In this paper, we first demonstrate that vdW-DF permits structural determination of molecular crystals and
second asses the character of the attractive dispersion
interactions. For the first, we calculate the lattice parameters for crystals of the cage-like hexamine, dodecahedrane, and cubane molecules. For the second, we compare an atom-centered asymptotic-1/r6 approximation,
frequently encountered in force-field methods6,30,31 and
in DFT-D, with the correlation energy provided by vdWDF. We also analyze consequences for DFT-D and question its performance at intermediate separations. The
particular choice of molecules was motivated in part by
their aesthetic appeal and in part for the property that
their symmetric geometries permits a simple assignment
and analysis of the strength of the asymptotic interactions even when expressed in terms of an atom-pair basis.
This paper has the following plan: The next section
presents the molecules hexamine, cubane, and dodecahedrane and their experimental crystal structures. The
third section deals with the computational details of
vdW-DF. The fourth gives the results for lattice parameters and bulk modulus. In the fifth section, we compare
the asymptotic 1/r6 -form of van der Waals interactions
with the full non-local correlation for different molecules,
and the vdW-DF potential energy curve with DFT-D cal-
2
FIG. 1: The molecular and crystal structure of hexamine, dodecahedrane, and cubane. The upper panels show the respective
molecular structures. The mid panels give schematics of the crystal structures, where the molecular symmetry and orientation
are highlighted by the use of a tetrahedron, dodecahedron, and a cube in place of the molecules. The lower panel shows, for
hexamine and dodecahedrane, facets of the crystal structure and for cubane the orientation of the cubane molecule in the unit
cell. Hexamine has a bcc unit cell, dodecahedron an fcc unit cell, and cubane a rhombohedral unit cell, with equal angles α
between the lattice vectors.
culations for a cubane dimer. The final section holds our
conclusions.
II.
PLATONIC MATTER
Plato asserted that the basic building blocks of Nature were five geometrical structures, today known as
platonic solids.32 These are the five convex polyhedra
with all faces, edges, and angles congruent. On the
molecular level, two of these solids have synthetic hydrocarbon (compounds of only carbon and hydrogen) analogues: cubane (C6 H6 ), which corresponds to the cube,
and dodecahedrane (C20 H20 ), which corresponds to the
twelve-faced dodecahedron. The fascination mathematicians held for these geometrical structures since antiquity
echoed in the more recent struggle to synthesize their hydrocarbon representatives. The highly strained bonds
posed the main challenge. Obstacles were eventually
overcome and in 1964 cubane was synthesized33 followed
by dodecahedrane in 1978.34 A third, tetrahedrane exists
only within a larger chemical structure35 and no crystallographic characterization exists. Tetrahedrane does
therefore not represent a good testing ground for vdWDF and we instead study the molecular crystal of hexamine (C6 H12 N4 ). Although not a hydrocarbon, it shares
several features with the platonic hydrocarbons: it is organic, its non-hydrogen atoms form a cage, and this cage
has the symmetry of platonic solid tetrahedron. The top
row of Fig. 1 shows the molecular structures (from left to
right) of hexamine, cubane and dodecahedrane, ordered
according to the ascendancy of their platonic analogues.
The second and third row of Fig. 1 show the crystal
structures. The first column shows hexamine, the mid,
cubane, and the final, dodecahedrane. The crystal structure of hexamine forms a bcc with a single molecule in
each unit cell and I4̄3m symmetry. Hexamine has been
used as a model system in numerous studies36–38 and the
crystal structure was determined as early as in 1923.39
The crystal structure of cubane is rhombohedral with
equal external angles and a single molecule in each cell.
The molecule is oriented according to the R3̄ space group;
3
a rotation about the [111] axis specifies the configuration
of the molecule.31 The crystal structure of dodecahedrane
is fcc, also with a single molecule in each cell and Fm3̄space group.40
Hexamine finds numerous industrial uses.41 For instance, it serves as a component in fuel tablets and as
an antibiotic.42 The platonic hydrocarbons find only hypothetical applications; cubane has been identified as a
potential high-energy fuel and explosive.43
In addition to being well-studied crystals, especially
hexamine and cubane, these crystals make attractive
testing grounds for vdW-DF for two more reasons. First,
their simple structures allow for brute-force determination of lattice parameters. Within the crystal symmetry,
this determination corresponds to mapping out the potential energy of a single parameter for hexamine and dodecahedrane, and three for cubane. The brute-force approach facilitates post-processing analysis, such as computation of bulk modulus,44 and makes it easier to evaluate how the choice of exchange in DFT affects crystal
structure and cohesion energy. Second, the high symmetry of the molecules reduces the large set of atom-to-atom
C6 coefficients to only a few equivalent values. This reduction simplifies the comparison of full vdW-DF calculations with atom-based asymptotic accounts of the van
der Waals forces.
III.
COMPUTATIONAL DETAILS
The crystal structure and bulk modulus of hexamine,
cubane, and dodecahedrane are determined with DFT
using a nonlocal density functional called vdW-DF. The
traditional semi-local GGA for exchange-correlation provides accurate bond lengths and charge density n(r),
but fails to capture correlated motion of separated electrons: the long-range dispersion forces. vdW-DF includes these correlations and can therefore account for
the structure and cohesion of sparse matter. Since details of the functional and its implementation are given
elsewhere,10,21,45–47 we focus mostly on computational
steps specific for determination of crystal structures.
The non-local correlation of vdW-DF takes the form
of a double-space integral,
Z
Z
1
nl
Ec [n] =
dr′ n(r) φ(r, r′ )n(r′ ) ,
(1)
dr
2 V0
V
over an interaction kernel φ(r, r′ ). Here V0 denotes the
central unit cell and V (formally) the entire space. The
kernel can be tabulated in terms of two parameters d and
d′ related to the local response q0 (r) and spatial separation |r−r′ | by d = q0 (r)|r−r′ | and d′ = q0 (r′ )|r−r′ |. The
remaining part of the exchange-correlation functional of
vdW-DF consists of the exchange part of revPBE48 and
the correlation of LDA:
vdWDF
Exc
= EcLDA + ExrevPBE + Ecnl [n] .
(2)
The total energy functional of vdW-DF, E vdWDF [n],
also includes the standard elecrostatic and kinetic-energy
terms within the Kohn-Sham scheme.49 It is convenient
to write this energy in terms of a semi-local part E0 [n]
containing all but the non-local correlation, so that
E vdWDF [n] = E0 [n] + Ecnl [n]. For input charge density
n(r), we use the result of semi-local calculations with
the PBE5 flavor of GGA. We will refer to calculations
with the PBE flavor of GGA as DFT-GGA. The charge
density could also have been obtained within fully selfconsistent vdW-DF.45,50 However, the two-step non-self
consistent procedure introduces only a slight approximation. Previous studies have documented that for systems
with small charge transfer, binding energies of non-self
consistent vdW-DF only differ by tiny amounts from fully
self-consistent energies.45,51
To speed up evaluation of the non-local correlation, we
introduce a radius cutoff based on the decay of van der
Waals forces at large separations. With this cutoff, the
kernel takes form
Z
Z
1
nl
Ec [n] ≈
dr′ φ(r, r′ )n(r′ ) . (3)
dr n(r)
2 V0
|r−r′ |<R
In the above expression, we see that the CPU-cost
for evaluating the non-local part of vdW-DF goes as
R3 O(V0 ). Thus, for large or periodic systems the computational costs increase linearly with system size. To cut
computational costs further, we introduce an extra radius
cutoff corresponding to the separation between dense and
sparse sampling of the charge-density grid. We note that
Ref. 50 reports a more elaborate scheme which considerably reduces CPU-costs, yet our simple measure was
sufficient for our non-self consistent calculations as the
underlying DFT-GGA calculations of electronic density
dominated time consumption.
The use of revPBE for exchange in vdW-DF was motivated by the fact that this exchange functional excludes unphysical binding effects at large distances.46,52
For a range of systems, vdW-DF overestimates binding
separations.21,46 Several studies indicate that this discrepancy can be attributed to the details of the exchange
functional.21,51,53–55 Puzder et al53 have demonstrated
that replacing Hartree-Fock exchange with revPBE improves binding separations for benzene dimers. Gulans
et al54 have shown that for a selected range of molecular
complexes the PBE exchange functional improves binding energies. We furthermore illustrate the sensitivity to
exchange by including results based on use of an alternative vdW-DF(PBE), where revPBE exchange has been
replaced by that of PBE. We do not argue for replacing
vdW-DF with vdW-DF(PBE), instead we simply explore
consequences of a different account of exchange.55
To calculate crystal parameters and cohesion energy,
we minimize the potential energy. In many respects this
vdW-DF structure determination is similar to those in
Refs. 2,23,47,56. The potential energy is given by the
difference between the total energy of the full crystal and
4
a reference energy for a system of isolated molecules
TABLE II: The vdW-DF prediction of lattice parameters, cohesion energy and bulk modulus for the crystals of hexamEcoh (a, {αθ}) = E
(a, {αθ})−E
(a → ∞, {αθ}) . ine, cubane and dodecahedrane compared with an alternative
(4)
vdW-DF(PBE) based on PBE-exchange and with experimental values. The experimental lattice parameters are based on
In the above equation the curly brackets are specific for
low temperature measurements, except for dodecahedrane.
the cubane crystal as it depends on three rather than one
vdWDF
vdWDF
parameter. In the reference calculation corresponding to
the reference energy, E vdWDF (a → ∞, {αθ}), the semilocal part is obtained somewhat differently from the part
containing the non-local correlation. For the former, we
effectively isolate the molecules in our periodic boundary
calculations, by using a unit cell of doubled size in all directions. This measure secures negligible charge overlap
between the molecules in the supercells. For the latter
(nonlocal) part we restrict the integral of Eq. 3 to the
central supercell to avoid coupling between the enlarged
unit cells. Hence, only non-local correlations within the
molecule contribute to the reference energy.
To enhance accuracy in the evaluation of the non-local
part of the potential energy, we systematically cancel a
small, but noticeable, grid dependence in the evaluation
of Eq. 3. This cancellation is performed by making sure
to use the same FFT grid spacing in the reference calculation as in the main calculation. Furthermore, we make
sure to place the isolated molecules in the same relative
configuration to the underlying grid in the reference calculation as in the main calculation. We thus perform an
additional reference calculation for every molecular configuration investigated. These measures have been used
to secure a high accuracy of the non-local part of vdWDF in several earlier studies.2,46,47,56
We map the potential energy landscape by varying the
lattice parameters of the molecular crystals within the
experimental crystal symmetry. The stiff cage-molecules
allow us to keep the internal coordinates of the molecules
frozen for all configurations. The molecular structures
are determined in isolation using the PBE flavor of GGA.
The resulting bond lengths will be compared with experimental data in the next section to verify the utility
of conventional DFT-GGA for the internal structure of
strained molecules.
TABLE I: Experimental and calculated bond lengths of hexamine, cubane, and dodecahedrane. The calculations were
done with the PBE flavor of GGA. l denotes the C-C bond
length for cubane and dodecahedrane and the C-N bond
length for hexamine. lCH denotes the carbon-hydrogen bond
length.
Parameter
l[Å]
lexp [Å]
lCH [Å]
exp
lCH
[Å]
a Ref.
58.
b Ref. 59.
c Ref 60.
Hexamine
1.472
1.476a
1.101
1.088a
Cubane
1.566
1.562b
1.095
1.097b
Dodecahedrane
1.549
1.544c
1.100
-
Parameter
vdW-DF
vdW-DF (PBE)
Exp.
a (Å)
Ecoh (eV)
B0 (GPa)
7.14
-1.0127
10.0
6.93
-1.427
14.0
6.910a
-0.827b
7.0c
a (Å)
α
θ
Ecoh (eV)
B0 (GPa)
5.45
73
47.5
-0.77
7.2
5.25
72.5
46.5
-1.15
14.8
5.20d
72.7d
46 d
-0.857e
-
a (Å)
Ecoh (eV )
B0 (GPa)
10.92
-1.46
12.2
10.56
-2.06
18.6
10.60f
Hex.
Cub.
Dod.
a Ref.
-
58. b Ref. 61. c Ref. 62. d Ref. 58. e Ref. 63. f Ref. 64.
The electronic-structure calculations rely on the planewave code DACAPO 57 using ultra-soft pseudopotentials.
In combination with a separate reference calculation as
previously discussed, we secure the convergence of the
non-local correlation by specifying an FFT-grid spacing
less than 0.13 Å. This spacing leads to an effective planewave energy cutoff of at least 500 eV. For all crystals, we
set the Monkhorst-Pack k-sampling to 4 × 4 × 4.
IV.
STRUCTURE DETERMINATION
A.
Molecular structure
Table I shows the calculated molecular structures.
It demonstrates DFT-GGA can account for the intramolecular bonding even in the highly strained cubane
molecules. The resulting bond lengths differ by less than
1 % from the experimental values. In first principle studies of molecular crystals, accurate determination of lattice parameters require accurate account of constituent
molecules.
Unlike many calculations, where determination of crystal structure often starts from the structure of the individual molecules, most experiments resolve the molecular structure by looking at the diffraction pattern of a
full molecular crystal. For our purposes of testing the
first-principle vdW-DF method, this is fortunate, as efforts to characterize molecular structures also generate
an abundance of experimental data on molecular crystal
structures.
5
FIG. 2: The potential energy curves for the crystal of hexamine and dodecahedrane (left panel) and corresponding contour plots
for the crystal of cubane (right panel). The left panel shows potential energy where the curves are normalized separately so that
the experimental lattice parameter and the calculated cohesion energy equals unity. The solid curve represents dodecahedrane
and the dashed hexamine. The difference in the optimal a/aexp value can be attributed to the experimental lattice parameter
being measured at low temperature for hexamine, but not for dodecahedrane. The right panel shows spline-interpolated contour
plots of the two-dimensional intersection of the three-dimensional potential energy landscape. The main figure corresponds
to the optimal value of the internal angle θ, while the insert corresponds to the optimal value of the unit cell length a. The
pronounced asymmetry of the curves and contours reflect the hard wall provided by Pauli repulsion.
B.
Crystal structure
Figure 2. shows the binding curves and contours found
by varying the molecular crystals’ lattice parameters
identified in Fig. 1. The simple crystal structures of
hexamine and dodecahedrane give a one-dimensional potential energy landscape, while the cubane crystal has a
three dimensional one. In the right panel of Fig. 2, we
display the αa and the αθ intersections, for the optimal
values of θ and a. The curves exhibit a pronounced asymmetry around their minimum. This asymmetry arise because the van-der Waals attraction is much softer than
the kinetic-energy repulsion.
Table. II contains the calculated results and experiment values. Standard vdW-DF performs well both for
lattice parameters and cohesion energy. If not directly
available, we obtain experimental cohesion energies by
correcting for gas phase and vibrational contributions to
the enthalpy of sublimation, using the method described
in Refs. 65,66.
The bulk modulus is obtained with use of polynomial
interpolation according to the scheme of Ziambaras and
Schröder.44 The required polynomial were constructed
using data from selected one-dimensional deformations.67
For hexamine, where the experimental bulk modulus is
available, the computed value show fair agreement with
the experimental value. There is also a trend for bigger
cage molecules to have a larger bulk modulus. We attribute this trend to the fact that for bigger molecules a
smaller relative part of the unit cell consists of soft intramolecular regions. Therefore, as the relative unit cell dimensions change, the distance between (stiff) molecules
changes more for big molecules than for small molecules.
For all three crystals, we find unit-cell volumes somewhat larger than the experimental ones. Similar overestimations have also been encountered in previous
studies.21,46,55 The alternative choice of PBE as exchange
functional, vdW-DF(PBE), influences results substantially. On one hand, it improves lattice parameters, almost to level of standard DFT for intra-molecular bonds.
On the other hand, the value of cohesive energy and bulk
modulus worsens. These results signal that the main discrepancy between experiments and vdW-DF stems from
the specific form of semi-local exchange.21,53–55
We note that an ab initio study of cubane has previously been performed at the LDA level.68 Being a local
functional, LDA has no physical basis for the van der
Waals binding that provides the cohesion of this molecular crystal. The spurious LDA binding arises from an
unphysical accounts of exchange.46,52 Once the molecular crystals are investigated with DFT-GGA, which has
an improved account of exchange, the binding essentially
vanishes.68
V.
ASYMPTOTIC PAIR POTENTIALS VERSUS
NON-LOCAL CORRELATION
The asymptotic van der Waals interactions between
two atoms or molecules goes as the simple power law
C6 /r6 , where C6 gives the strength of the interaction.
This familiar result can be derived from second-order perturbation theory,69 or from an analysis of the shifts in the
zero-point motions of the electron.70 A common strategy
in force-field methods and empirical extensions of DFT
6
than to the surface atoms. Within vdW-DF, Kleis et al23
demonstrated that for interactions in nanotube bundles,
the force stems primarily from the electron tail around
the nanotube, and that as the tubes get closer, higher order moments dominate over the asymptotic interaction.
Part of this enhancement can be interpreted as an imageplane effect.
Here we investigate whether an atom-centered 1/r6
form is a good approximation for the non-local correlation of molecular dimers. As our argument is based
on the asymptotic vdW-DF account of van der Waals
forces, we also discuss for these molecules the accuracy
of the asymptotic account. In the first subsection, we
compare the non-local correlation of vdW-DF for dimers
of hexamine, cubane, and dodecahedrane with its corresponding atom-centered pair potentials. We document
a significant enhancement of the non-local correlation at
short (binding) separations and at intermediate separations, one to three Ångström beyond typical binding separations. In the second subsection, we discuss the accuracy of our C6 coefficients. In the third, we investigate
consequences of this result, in particular for the use of
DFT-D. We will argue that although standard DFT-D
methods can provide good descriptions of both short and
asymptotic separations, their asymptotic atom-centered
form does not describe the enhancement of correlation at
intermediate separations exhibited by vdW-DF.
FIG. 3: Comparison between the non-local correlation of
vdW-DF and atomic pairs potentials generated with asymptotic vdW-DF for a dimer of hexamine, dodecahedrane, and
cubane (from top to bottom) in configurations corresponding to nearest neighbors in their respective crystal. The
dashed curve gives the APP-vdWDF result, while the solid
curve gives non-local correlation of vdW-DF: ∆Ecnl (d) =
Ecnl (d) − Ecnl (d → ∞). The horizontal axis gives the difference from the vdW-DF crystal binding separation d0 , which
respectively takes the values 6.18, 5.45 and 7.72 Å. The difference between the two curves demonstrates the enhancement
of non-local correlations over the vdW-DF asymptotic atombased account at relevant binding separations and intermediate separations.
A.
Asymptotic vdW-DF
For the the asymptotic van der Waals forces, the C6
coefficient between two fragments, A and B, can be computed from the general formula
Z
3 ∞
C6AB =
du αA (iu)αB (iu) ,
(5)
π 0
where α(ω) is the polarizability of the fragment. In order
to calculate the C6 coefficients we approximate α(ω) with
the local external-field susceptibility of vdW-DF,
χvdW−DF
(ω, r) =
A
is to adopt such an asymptotic form at all separations in
terms of atom-centered pair-potentials (APP). However,
we can not take for granted that the asymptotic behavior
should hold for separations closer to that of intramolecular binding. On the contrary, because van der Waals
forces arise from correlated motion of electrons and not
from the atomic nuclei, several mechanisms affect and
enhance the interaction at short and intermediate separations: higher order moments contribute, polarizability
changes as charge is distorted and finite-k dispersion of
the electronic response becomes important. Zaremba and
Kohn71 considered adsorption of noble atoms on surfaces
and documented a significant enhancement of dispersion
energies over an atom-centered account; their asymptotic
1/d3 form use the distance to an image plane d, rather
nA (r)
2
[9q0 (r)2 /8π] − ω 2
,
and a polarizability given by
Z
vdW−DF
αA
(ω) = d3 r χvdW−DF
(ω, r) .
A
(6)
(7)
To generate vdW-DF based atom-centered pair potentials (APP-vdWDF),
P we first partition the full charge
density n(r) =
i ni (r) among the atoms of the
molecules with aid of Bader analysis.72,73 Based on this
charge partition, we calculate the atom-to-atom C6 coefficients using Eq. (5). Initially, this procedure generates
N 2 different C6 coefficients for a molecular dimer of N
atoms per molecule. For a dodecahedrane dimer, we get
as much as 40 × 40 = 1600 coefficients. Fortunately, because of the high symmetry of the isolated molecules,
7
this number reduces to only three equivalent values for
cubane and dodecahedrane: C6C−C , C6C−H , and C6H−H .
For hexamine, the extra nitrogen atoms lead to six coefficients. As noise in the electronic density affects the
value of the coefficients, we average over a large set of
equivalent values to obtain the final values.
TABLE III: Computed values of the C6 coefficients (Hartree
atomic units) for different pairs of atoms within respective
molecules, using the asympotic form of vdW-DF with charge
density as partitioned with a Bader analysis. The figure also
shows the molecule-molecule C6 coefficients obtained with
the Andersson-Langreth-Lundqvist (ALL) scheme,74 and the
molecule-molecule C6 coefficents obtained with the use of parameters given in DFT-D schemes.16–18
C6vdWDF
C-C
C-H
H-H
N-N
N-C
N-H
mol-mol
mol-mol
mol-mol
mol-mol
mol-mol
(vdW-DF)
(ALL)
(Wu)
(Grimme)
(Jurečka)
hex
4.44
4.04
3.98
32.6
12.0
11.0
3470
3270
4340
4000
4790
cub
13.5
6.87
3.72
1990
1940
2600
2630
3120
dod
11.2
6.45
4.07
11300
10200
16300
16400
19500
The upper part of Table III shows the C6 coefficients
calculated within asymptotic vdW-DF, both as partitioned according to the Bader analysis and as evaluated
for the entire molecule. The calculated coefficients per
atomic pair deviates much from a naive assignment of
the full coefficient of the molecule according to the number of valence electrons of the underlying atom. In such a
scheme the C-C coefficient would be 16 times larger than
the H-H coefficient, while in fact, for cubane and dodecahedrane, it is merely three-four times stronger. This result can be attributed to the relative stronger response of
the low-density regions surrounding the hydrogen atoms
as q0 (r) ∝ [n(r)]1/3 (in the homogeneous limit), and these
areas dominate Ecnl . The somewhat anomalous values for
hexamine can be attributed to our Bader analysis scheme
partitioning a significant portion of the charge density
near the carbon atoms to the centrally located nitrogen
atoms.75 For cubane and dodecahedrane the partitioning
was similar, and the differences in atomistic C6 values
show that they are influenced by their local enviroment.
Having generating C6 coefficients appropriate for a
comparison between the asymptotic account and the
full correlation of vdW-DF, we study dimers of hexamine, cubane and dodecahedrane at different separations.
We choose orientations that are given by the nearestneighbor configurations in the crystals. The total asymp-
totic non-local correlation of APP-vdWDF reads
vdWDF
=
Eapp
XX
i
j
C6ij
,
|ri − r′j |
(8)
where i and j label atoms in separate molecules of the
dimer.
Figure 3 shows results for the full non-local correlation of vdW-DF (solid curve) and that of APP-vdWDF
(dashed curve). At large separations these two curves
converge. In contrast, they differ significantly at relevant
binding (short) separations and at intermediate separations. At these separations the non-local correlation is
almost double as large as that of APP-vdWDF, for both
hexamine and cubane, while for the biggest molecule, dodecahedrane, the difference is smaller. For all dimers, we
also need to go to relatively large separations to recover
asymptotic values for the non-local correlation.
As vdW-DF is an ab initio functional, based on a set of
exact sum rules,10 the strong enhancement of non-local
correlations at short and intermediate separations indicates that the asymptotic form neglects important contributions. It also highlights that in constructing modified
semi-empirical van der Waals functionals,11,12 a fitting
of the asymptotic functional to C6 coefficients does not
guarantee an accurate description at short binding separations.
B.
Comparison of C6 coefficients
The lower part of of Table III, gives the molecular C6 coefficients calculated with asymptotic vdW-DF,
the Andersson-Langreth-Lundqvist (ALL) scheme,74 and
computed using the atomic coefficients of Wu et al,16
Grimme,17 and Jurečka et al.18 The ALL and vdWDF give coefficients smaller than that used in DFT-D
schemes.
We expect that Grimme and Wu provide good values
for molecular C6 coefficients, because their underlying
atomistic coefficients were fitted to reproduce a range
of accurate molecular C6 coefficients calculated from experimental molecular polarizabilities (Ref. 16 and references therein). The coefficients of Jurečka,18 give somewhat larger molecular C6 coefficients. This comes from
the use of Slater-Kirkwood average76 for C6 coefficients
between different atomic species, while keeping those of
Grimme for indentical atomi species (the C6 coefficients
of Grimme are optimized for a different average.
Asymptotic vdW-DF and ALL likely underestimates
the C6 coefficients for these molecules; they are 20-30 %
smaller than that used in DFT-D methods. For the similar ALL scheme, Ref. 77 reports an underestimation of
C6 coefficients for larger molecules, in particular for benzene and C60. A difference between the ALL scheme and
asymptotic vdW-DF is that for the former a hard cutoff
accounts for plasmon damping, while for the latter, the
local response q0 (r) provides a smooth cutoff. There is
8
good consistency between the two methods. Both methods also assume a local, scalar, relationship between the
applied and the full electric field, which is an approximation for finite-sized objects.74 We speculate that this approximation contributes to the underestimation of C6 coefficient for the investigated, relatively large, molecules.
TABLE IV: Atomistic C6 coefficients for cubane used in APPvdWDF (calculated with asymptotic vdW-DF and charge
density partitioned with Bader analysis), and coefficients used
in DFT-D methods.
C6vdWDF
C-C
C-H
H-H
a Ref.
b Ref.
c Ref.
vdW-DF
13.5
6.8
3.7
Wua
22.06
7.89
2.83
. . .
Grimmeb
28.3
5.01
2.75
Jurečkac
28.3
8.82
2.75
16
17
18
Table IV shows atomistic C6 coefficients for cubane as
calculated with asymptotic vdW-DF, and given by Wu,16
Grimme,17 and Jurečka,18 for use in DFT-D schemes.
vdW-DF weights the relative response of the carbon less
than that of the hydrogen, compared with the values
of DFT-D.75 This property could relate to the above
mentioned approximate treatment of electrodynamics. It
could also relate to the carbon atoms being located somewhat inside the molecule, having a different local charge
density and responding less to external fields than an
atom on the exterior would; in contrast, DFT-D does
not discriminate between atoms at different locations.
C.
Consequences for atom-based pair potentials:
the binding curve of cubane
The vdW-DF results presented in the first subsection
shows that a simple asymptotic account only partially
captures correlation effects at short and intermediate separations. As DFT-D use such an asymptotic form to
describe non-local correlations, this result stands in apparent contrast to the many successful applications of
DFT-D.16–18,78–80
To understand consequences of our result for methods such as DFT-D, we must first consider other effects
that could contribute to the difference between APPvdWDF and the non-local correlation of vdW-DF. In
vdW-DF, correlations are described by EcLDA + Ecnl and
hence Ecnl also accounts for semi-local correlations.10,46
Second, we must consider the specific designs of actual
DFT-D schemes, because these could counteract the lack
of enhancement of non-local correlations. To describe
exchange-correlation, DFT-D combines the asymptotic
atom-centered form with a semi-local GGA account. It
also introduces fitting parameters to be used in combination with a specific GGA flavors.
The top panel of Fig. 4 shows that gradient corrections
does not account for the difference between the full vdWDF and the APP-vdWDF results. For APP-vdWDF,
we can combine the purely semi-local correlation of PBE
with APP-vdWDF (in a new description APP-mod),
vdWDF
vdWDF
EAPP−mod
(d) = EAPP
(d) + ∆EcPBE (d) − ∆EcLDA (d) ,
(9)
to assess the magnitude of purely semi-local corrections
relative to the difference between vdW-DF and APPvdWDF. In this APP account, we have introduced the
LDA and PBE terms: ∆Ec (d) = Ec (d) − Ec (d → ∞).
We focus our discussion on the cubane dimer. The
lower thin solid curve gives APP-mod; the thick solid
curve gives the non-local correlation of vdW-DF. The
curves are shown as a function of the intermolecular distance d, with d0 indicating the vdW-DF binding distance
in the crystal. At short separations, the thin curve lies
closer to the thick curve than the corresponding APPvdWDF result (thick dashed curve). Thus, some of the
difference between APP-vdWDF and the non-local correlation arises from a lack of semi-local correlation contributions in APP-vdWDF.10 However, even with this
inclusion the difference is still significant, and at intermediate separations it remains undiminished.
The middle panel of Fig. 4 details the effects of using
the semi-empirical fitting of DFT-D and the larger C6
coefficients. We compare the binding curve for cubane
obtained by use of vdW-DF with the binding curve obtained with DFT-D calculations. We select the schemes
of Grimme17 and Jurečka,18 as these provide parameters for the PBE flavor of exchange-correlation, which is
available to us.81
The DFT-D scheme of Grimme scales the strength
of the dispersive interaction to a particular semi-local
exchange-correlation to achieve good performance at
short separations. For the PBE flavor of GGA, the cohesion energy of a dimer reads
DFT−D,G
PBE
Ecoh
= Ecoh
+sPBE
6
X
ij
fPBE,G
(|ri −rj |)
ij
C6ij
,
|ri − rj |6
(10)
where s6 = 0.7. This scheme therefore sacrifices the
asymptotic description in favor of the a good description
of binding separations. The scheme of Jurečka instead
adjust parameters of the damping function to the flavor
of semi-local exchange-correlation,
DFT−D,J
PBE
Ecoh
= Ecoh
+
X
ij
ij
fPBE,J
(|ri − rj |)
C6ij
,
|ri − rj |6
(11)
and ensure that
→ 1 for large separations. The
DFT-D scheme can therefore, in principle, describe both
asymptotic and binding separations. However, the question remains on how it performs for systems where characteristic separations lies between these two limits.
The binding curves of vdW-DF(revPBE) (the solid
curve) and for the DFT-D scheme of Grimme(PBE) (upij
fPBE,J
9
responding DFT-D scheme binds stronger than vdW-DF.
The lower panel of Fig. 4 shows that the difference at
intermediate separations comes primarily from the varying accounts of correlation. The solid curve shows the
non-local correlation of vdW-DF. The dash-dotted upper (lower) curve shows
DFT−D,G(J)
Eapp
=
X
ij
fPBE,G(J)
(|ri − rj |)
ij
C6ij
|ri − rj |6
+ ∆EcPBE (d) − ∆EcLDA (d) ,
FIG. 4: Comparison between different accounts of non-local
correlation for a dimer of cubane. In the upper panel, the
dashed curve gives APP-vdWDF. The thin solid curve gives
the sum of APP-vdWDF and the correlation of PBE. The
thick solid curve gives the non-local correlation of vdW-DF.
In the middle panel, the thick curve gives the cohesion envdWDF
ergy using vdW-DF, Ecoh
(d), the upper (lower) dashdotted curve gives the DFT-D binding curve as given by
DFT−D,G
DFT−D,J
Grimme, Ecoh
(d) (Jurečka, Ecoh
(d)). The lower
panel gives the corresponding non-local correlation of vdWDF and Grimme (Jurečka). The thin solid curve gives the difference between exchange of revPBE and PBE. The horizontal axis gives the difference from the vdW-DF crystal binding
separation d0 for cubane.
per dash-dotted curve) and Jurečka(PBE) (lower dashdotted curve) indicates an underestimation of DFT-D
at intermediate separations. The minimum at negative
d − d0 shows that use of DFT-D would improve lattice
constants over vdW-DF.55 The two DFT-D schemes give
quite differing binding energies. Both results show that
a binding energy at the same magnitude as vdW-DF can
be achieved, even with an atom-based asymptotic form
of the attractive potential. The energy of vdW-DF is significantly larger at intermediate separations than that of
Jurečka despite that C6 coefficients of asymptotic vdWDF are underestimated (while those of Jurečka are likely
to be somewhat overestimated) and despite that the cor-
(12)
which for DFT-D corresponds best to the non-local correlations provided by ∆Ecnl in vdW-DF. The figure also
shows that the difference partly cancels, at short but
not at intermediate separations, with the energy difference between the exchange flavors of revPBE and PBE,
δEx = ∆ExrevPBE − ExPBE . Thus, for certain exchange
functionals, adding an asymptotic atom-based account of
non-local correlation can generate good binding values,
yet our vdW-DF results indicate that this framework is
not optimal for describing interactions at intermediate
separations.
Our results suggest that APPs and DFT-Ds could be
improved at intermediate separations. A possible strategy is to replace the atomic separations r by an effective
separation r−r0 , where r0 reflects the image-planes found
for nanotubes and surfaces in Refs. 23,71. Keeping this
(surface-physics) effect would increase the strength of the
dispersion interactions at shorter separations.
In summary, the results of this section suggests that
an asymptotic atom-based pair potentials has a limited
transferability over the full range of separations. Thus,
for schemes using such a form, our results raises questions on their ability to generate accurate results under
broad condition (having multiple characteristic separations), for instance involving phase transitions or processes that drive the system out off equilibrium, in protein unfolding, in phase transitions, or simply for systems which have competing interactions.1 We note that
there is no guarantee that vdW-DF, in its current form,
can provide an accurate account under such broad conditions. A vdW-DF limitation is here exemplified by the
likely underestimation of C6 coefficients for the platonic
molecules. Nevertheless, we argue that non-local functionals, like vdW-DF, hold the most promise for dealing with molecular configurations under broad conditions. This is because an electron-based approach provides a framework which naturally includes image-plane
and multipole effects. It therefore holds the key to an account which describe the variation in dispersive response
over the full range of separations.
VI.
CONCLUSIONS
For the three cage-like organic molecular crystals, hexamine, cubane and dodecahedrane, vdW-DF gives lattice
10
parameters and cohesion energy that agree well with experiments, although for all three crystals the unit cell
volumes are overestimated. A substantial sensitivity
of lattice parameters and cohesion energy to the flavor
of semi-local exchange signals that this overestimation
stems mostly from the chosen form of exchange functional.
We have also shown that, at short and intermediate
separations, the full non-local correlation of vdW-DF is
considerably larger than its corresponding atom-based
asymptotic account. Notwithstanding that the asymptotic account of vdW-DF likely needs improvement, this
enhancement indicates that the asymptotic 1/r6 form of
atomic pair potentials, by construction, can not give a
transferable account over a large range of separations.
This paper underlines the usefulness of studying simple model sparse systems as molecular crystals to gain
insight into methods intended for the study of sparse
and supramolecular systems.1 Both DFT-D and vdWDF benefit from such testing, because they are designed
to be parameter-free, and an accurate account of molec-
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
J.M. Lehn, Science 260, 1762 (1993).
J. Kleis, B. I. Lundqvist, D. C. Langreth, E. Schröder,
Phys. Rev. B 76, 100201 (2007).
S. H. Vosko, L. Wilk and M. Nusair, Can. J. Phys. 58,
1200 (1980).
J. P. Perdew and Y. Wang, Phys. Rev. B 45, 13244 (1992).
J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev.
Lett. 77, 3865 (1996).
Day, G. M., Chisholm, J., Shan, N., Motherwell, W. D. S.
and Jones, W, Cryst. Growth Des.4, 1327 (2004).
W. Kohn, Y. Meir, and D. E. Makarov, Phys. Rev. Lett.
80, 4153 (1998).
J.F. Dobson, and J. Wang, Phys. Rev. Lett. 82, 2123
(1999).
H. Rydberg, M. Dion, N. Jacobson, E. Schröder, P.
Hyldgaard, S. I. Simak, D. C. Langreth, and B. I.
Lundqvist, Phys. Rev. Lett. 91, 126402 (2003).
M. Dion, H. Rydberg, E. Schröder, D. C. Langreth, and B.
I. Lundqvist, Phys. Rev. Lett. 92 (2004); Phys. Rev. Lett.
95, 109902(E) (2005).
O. A. Vydrov and T. Van Voorhis, J. Chem. Phys. 130,
104105 (2009).
O. A. Vydrov and T. Van Voorhis, Phys. Rev. Lett. 103,
063004 (2009).
M. A. Basanta, Y. J. Dappe, J. Ortega, and F. Flores,
Europhys. Lett. 70, 355 (2005).
Y. J. Dappe, M. A. Basanta, F. Flores, and J. Ortega,
Phys. Rev. B 74, 205434 (2006).
X. Wu, M. C. Varga, S. Nayak, V. Lotrich, G. Scoles, J.
Chem. Phys. 115, 8748 (2001).
Q. Wu and W. Yang, J. Chem. Phys. 116, 515 (2002).
S. Grimme, J. Comput. Chem. 25, 1463 (2004).
P. Jurečka, J. Černy, P. Hobza, D. R. Salahub, J. Comput.
Chem. 28, 555 (2007).
V. Barone, M. Casarin, D. Forrer, M. Pavone, M. Sambi,
A. Vittadini, J. Comp. Chem. 30, 934 (2008).
ular crystals would suggest an accurate account of more
complex assemblies of similar molecules. The presented
molecular crystals provide particularly accessible cases of
extended systems and can therefore be used in conjunction with future development of pair-potential methods
and exchange-correlation functionals.
VII.
ACKNOWLEDGEMENTS
We acknowledge E. Schröder for useful comments and
discussions. We thank I. Sinno for help in designing the
molecular-crystal schematics and Ø. Borck for access to
his code for Bader analysis. The Swedish National Infrastructure for Computing (SNIC) is acknowledged for
computer allocation and for KBs participation in the
graduate school NGSSC. The work was supported by the
Swedish research Council (Vetenskapsrådet VR) under
621-2008-4346.
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
M. A. Neumann and M. Perrin. J. Phys. Chem. B 109,
15531 (2005).
D. C. Langreth, B. I. Lundqvist, S. D. Chakarova-Kack, V.
R. Cooper, M. Dion, P. Hyldgaard, A. Kelkkanen, J. Kleis,
L. Kong, S. Li, P. G. Moses, E. Murray, A. Puzder, H. Rydberg, E. Schrder, and T. Thonhauser, J. Phys.: Condens.
Matter 21, 084203 (2009).
V.R. Cooper, T. Thonhauser, A. Puzder, E. Schröder, B.I.
Lundqvist, and D.C. Langreth, J. Am. Chem. Soc. 130,
1304 (2008).
J. Kleis, E. Schröder, and P. Hyldgaard, Phys. Rev. B 77
205422 (2008).
A. Kelkkanen, B. I. Lundqvist, and J. K. Nørskov, J.
Chem. Phys. 131, 046102 (2009).
L. Kong, G. Román-Pérez, J. M. Soler, and D. C. Langreth, Phys. Rev. Lett. 103, 096103 (2009).
P. Sony, P. Puschnig, D. Nabok, and C. Ambrosch-Draxl,
Phys. Rev. Lett. 99, 176401 (2007).
N. Atodiresei, V. Caciuc, P. Lazic, S. Blügel, Phys. Rev.
Lett. 102, 136809 (2009).
K. Berland, T. L. Einstein, and P. Hyldgaard, Phys. Rev.
B 80, 155431 (2009).
K. Toyodaa, Y. Nakanoa, I. Hamadaa, K. Leec, S. Yanagisawaa, and Y. Morikawaa, Surf. Sci. 603, 2912 (2009).
A. Gavezzotti, Modell. Simul. Mater. Sci. Eng. 10, R1
(2002).
T. Yildirim, P. M. Gehring, D. A. Neumann, P. E. Eaton,
and T. Emrick, Phys. Rev. Lett. 78, 4938 (1997).
Plato, Timaeus, http://classics.mit.edu/Plato/timaeus.html.
P. E. Eaton, and T. W. Cole, J. Am. Chem. Soc., 86, 3157
(1964).
L. A. Paquette, R. J. Ternansky, D. W Balogh, G. J. Kentgen, J. Am. Chem. Soc. 105, 5446 (1983).
M. N. Glukhovtsev, S. Laiter, A. Pross, J. Phys. Chem.
99, 6828 (1995)
M. S. Ahmed, Acta Cryst. 5, 587 (1952).
11
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
R. Rao, T. Sakuntala, S.K. Deb, A.P. Roy, V. Vijayakumar, B.K. Godwal, and S.K. Sikka, Chem. Phys. Lett. 313,
749 (1999).
G. H. Heilmeier, App. Opt. 3, 1281 (1964).
R. G. Dickison and A. L. Raymond, J. Am. Chem. Soc.
45, 22 (1923).
B. S. Hudson, D. G. Allis, S. F. Parker, A. J. RamirezCuesta, H. Herman, and H. Prinzbach, J. Phys. Chem. A
109, 3418 (2005).
A. Almdari and F. Tabkhi, Chem. Eng. and Proc. 43, 803,
(2003).
Ullmann’s Encyclopedia of Industrial Chemistry 2008.
S. Borman, Chem. Sci. News 72, 34 (1994).
E. Ziambaras and E. Schröder, Phys. Rev. B 68, 064112
(2003).
T. Thonhauser, V. R. Cooper, S. Li, A. Puzder, P.
Hyldgaard, and D. C. Langreth, Phys. Rev. B 76, 125112
(2007).
D. C. Langreth, M. Dion, H. Rydberg, E. Schröder, P.
Hyldgaard, B. I. Lundqvist, Int. J. Quantum Chem. 101,
599 (2005).
E. Ziambaras, J. Kleis, E. Schröder, and P. Hyldgaard,
Phys. Rev. B 76, 155425 (2007).
Y. Zhang and W. Yang, Phys. Rev. Lett. 80, 890 (1998).
W. Kohn and L. Sham, Phys. Rev. 140, A1133 (1965).
G. Román-Párez, and J. M. Soler Phys. Rev. Lett. 103,
096102 (2009).
S. Li, V. R. Cooper, T. Thonhauser, A. Puzder, and D.C.
Langreth, J. Phys. Chem. A 112, 9031 (2008).
X. Wu, M. C. Vargas, S. Nayak, V. Lotrich, and G. Scoles,
J. Chem. Phys. 115, 8748 (2001).
A. Puzder, M.Dion, and D. C. Langreth, J. Chem. Phys.
124, 164105 (2006).
A. Gulans, M. J. Puska, and R. M. Nieminen, Phys. Rev.
B 79, 201105 (2009).
Several new exchange flavors designed for use within the
vdW-DF framework promises to improve the vdW-DF
account of binding separations: V. R. Cooper, arXiv:
0910.1250v1; J. Klimeš, D. R. Bowler, and A. Michaelides
J. Phys.: Condens. Matter 22, 022201 (2010).
S. D. Chakarova-Käck, E. Schröder, B. I. Lundqvist, and
D. C. Langreth, Phys. Rev. Lett. 96, 146107 (2006).
Opensource
code
DACAPO,
http://www.fysik.dtu.dk/CAMPOS/
L.N. Becka and D.W.J. Cruickshank, Proc. Roy. Soc. A
273 (1963).
L. Hedberg, K. Hedberg, P.E. Eaton, N. Nodari and A.G.
Robiette. J. Am. Chem. Soc. 113, 1514 (1991).
K. Schmidt, M. Springborg, M. Bertau, F. Wahl, A.
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
Weiler, K. Scheumann, J. Worth, M. Keller M, H.
Prinzbach, Tetrahedron 53, 10029 (1997).
E. A. Arnautova, M. V. Zakharova, T. S. Pivina, E. A.
Smolenskii, D. V. Sukhachev, and V. V. Shcherbukhin,
Rus. Chem. Bul. 45, 1066 (1996).
G. N. Ramachandran and W. A. Wooster, Acta Cryst. 4,
431 (1951).
J. S. Chickos and W. E. Acree, J. Phys. Chem. Ref. Data
31, 537 (2002).
J. C. Galluci, C. W. Doecke, and L. A. Paquette, J. Chem.
Phys. 108, 1343 (1986).
D. Nabok, P. Puschnig, and C. Ambrosch-Draxl, Phys.
Rev. B 77, 245316 (2008).
A. T. Hagler, E. Huler, and S. Lifson, J. Am. Chem. Soc.
96, 5319 (1974).
W.F. Perger, J. Criswell, B. Civalleri, and R. Dovesi,
Comp. Phys. Comm. 180, 1753 (2009).
T. Yildirim, S. Ciraci, Ç. Kiliç, and A. Buldum, Phys. Rev.
B 62, 7625 (2000).
H. Margenau and N. R. Kestner, Theory of Interatomic
Forces (Pergamon Press, Oxford, 1969), 2nd ed.
G.Mahan, J. Chem. Phys. 43, 1569 (1965).
E. Zaremba and W. Kohn, Phys. Rev. B 13, 2270 (1976).
R. F. W. Bader, Atoms in Molecules: A Quantum Theory
(Oxford University Press, Oxford, 1990).
G. Henkelman, A. Arnaldsson, and H. Jónsson, Comput.
Mater. Sci. 36, 254 (2006).
Y. Andersson, D.C. Langreth, and B. I. Lundqvist, Phys.
Rev. Lett. 76, 102 (1996).
We have explicitly verified that the precise assignment of
the total molecular response onto the atoms does not alter the APP-vdWDF curves significantly. An assignment
of the total asymptotic vdW-DF response according to the
relative strength of Grimme yield a weak increase in the
difference between APP-vdWDF and the non-local correlation.
J. C. Slater and J. G. Kirkwood, Phys. Rev. 37, 682 (1931).
Y. Andersson and H. Rydberg, Phys. Scr. 60, 211 (1999).
M. Pavone, N. Rega and V. Barone, Chem. Phys. Lett.
452, 333 (2008).
Claudio Morgado, Mark A. Vincent, Ian H. Hillier and
Xiao Shan, Phys. Chem. Chem. Phys. 9, 448 (2007).
S. Grimme, J. Antony, T. Schwabe and C. MckLichtenfeld, Org. Biomol. Chem. 5, 741, (2007).
In our comparison, we neglect that these schemes are optimized for a particular basis set, arguing that our comparison does not require a fine-tuned version of DFT-D.