A major development in DFT during the mid-2000s was the recognition that, first
of all, semi-local density functionals do not properly capture dispersion (van
der Waals) interactions, a problem that has been addressed only much more
recently by the non-local correlation functionals discussed in
Section 5.7.1; and second, that a cheap and simple solution to this
problem is to incorporate empirical potentials of the form -C6/R6, where
the C6 coefficients are pairwise atomic parameters. This approach, which is
an alternative to the use of a non-local correlation functional, is known as
dispersion-corrected DFT (DFT-D).
459
Wiley Interdiscip. Rev.: Comput. Mol. Sci.
(2011),
1,
pp. 211.
Link
,
453
Chem. Rev.
(2016),
116,
pp. 5105.
Link
There are currently three unique DFT-D methods in Q-Chem. These are requested via the $rem variable DFT_D and are discussed below.
DFT_D
DFT_D
Controls the empirical dispersion correction to be added to a DFT calculation.
TYPE:
LOGICAL
DEFAULT:
None
OPTIONS:
FALSE
(or 0) Do not apply the DFT-D2, DFT-CHG, or DFT-D3 scheme
EMPIRICAL_GRIMME
DFT-D2 dispersion correction from Grimme
458
J. Comput. Chem.
(2006),
27,
pp. 1787.
Link
EMPIRICAL_CHG
DFT-CHG dispersion correction from Chai and Head-Gordon
207
Phys. Chem. Chem. Phys.
(2008),
10,
pp. 6615.
Link
EMPIRICAL_GRIMME3
DFT-D3(0) dispersion correction from Grimme (deprecated as
of Q-Chem 5.0)
D3_ZERO
DFT-D3(0) dispersion correction from Grimme et al.
450
J. Chem. Phys.
(2010),
132,
pp. 154104.
Link
D3_BJ
DFT-D3(BJ) dispersion correction from Grimme et al.
452
J. Comput. Chem.
(2011),
32,
pp. 1456.
Link
D3_CSO
DFT-D3(CSO) dispersion correction from Schröder et al.
1136
J. Chem. Theory Comput.
(2015),
11,
pp. 3163.
Link
D3_ZEROM
DFT-D3M(0) dispersion correction from Smith et al.
1190
J. Phys. Chem. Lett.
(2016),
7,
pp. 2197.
Link
D3_BJM
DFT-D3M(BJ) dispersion correction from Smith et al.
1190
J. Phys. Chem. Lett.
(2016),
7,
pp. 2197.
Link
D3_OP
DFT-D3(op) dispersion correction from Witte et al.
1378
J. Chem. Theory Comput.
(2017),
13,
pp. 2043.
Link
D3
Automatically select the “best” available D3 dispersion correction
D4
DFT-D4 dispersion correction from Caldeweyher et al.
162
J. Chem. Phys.
(2017),
147,
pp. 034112.
Link
,
163
J. Chem. Phys.
(2019),
150,
pp. 154122.
Link
,
164
Phys. Chem. Chem. Phys.
(2020),
22,
pp. 8499.
Link
RECOMMENDATION:
Use D4 if the specified functional is avialable. Currently, only a subset of functionals in DFT-D4 is supported.
It includes B3LYP, B97, B1LYP, PBE0, PW6B95, M06L, M06, WB97, WB97X, CAMB3LYP, PBE02, PBE0DH, MPW1K, MPWB1K, B1B95, B1PW91, B2GPPLYP, B2PLYP, B3P86, B3PW91, O3LYP, REVPBE,
REVPBE0, REVTPSS, REVTPSSH, SCAN, TPSS0, TPSSH, X3LYP, TPSS, BP86, BLYP, BPBE, MPW1PW91, MPW1LYP, PBE, RPBE, and PW91.
The oldest of these approaches is DFT-D2,
458
J. Comput. Chem.
(2006),
27,
pp. 1787.
Link
in which the
empirical dispersion potential has the aforementioned form, namely, pairwise
atomic -C/R6 terms:
ED2disp=-s6atoms∑Aatoms∑B<A(C6,ABR6AB)fD2damp(RAB). | (5.23) |
This function is damped at short range, where R-6AB diverges, via
fD2damp(RAB)=[1+e-d(RAB/R0,AB-1)]-1 | (5.24) |
which also helps to avoid double-counting of electron correlation effects, since short- to medium-range correlation is included via the density functional. (The quantity R0,AB is the sum of the van der Waals radii for atoms A and B, and d is an additional parameter.) The primary parameters in Eq. (5.23) are atomic coefficients C6,A, from which the pairwise parameters in Eq. (5.23) are obtained as geometric means, as is common in classical force fields:
C6,AB=(C6,AC6,B)1/2 | (5.25) |
The total energy in DFT-D2 is of course EDFT-D2=EKS-DFT+ED2disp.
DFT-D2 is available in Q-Chem including analytic gradients and frequencies,
thanks to the efforts of David Sherrill’s group. The D2 correction can be used
with any density functional that is available in Q-Chem, although its use
with the non-local correlation functionals discussed in Section 5.7.1
seems inconsistent and is not recommended. The global parameter s6 in
Eq. (5.23) was optimized by Grimme for four different
functionals,
458
J. Comput. Chem.
(2006),
27,
pp. 1787.
Link
and Q-Chem uses these as the default values:
s6=0.75 for PBE, s6=1.2 for BLYP, s6=1.05 for BP86, and s6=1.05 for B3LYP. For all other functionals, s6=1 by default. The D2
parameters, including the C6,A coefficients and the atomic van der Waals
radii, can be modified using a $empirical_dispersion input section. For example:
$empirical_dispersion S6 1.1 D 10.0 C6 Ar 4.60 Ne 0.60 VDW_RADII Ar 1.60 Ne 1.20 $end
Values not specified explicitly default to the values optimized by Grimme.
Note: 1. DFT-D2 is only defined for elements up to Xe. 2. B97-D is an exchange-correlation functional that automatically employs the DFT-D2 dispersion correction when used via METHOD = B97-D.
An alternative to Grimme’s DFT-D2 is the empirical dispersion correction of
Chai and Head-Gordon,
207
Phys. Chem. Chem. Phys.
(2008),
10,
pp. 6615.
Link
which uses the same form as
Eq. (5.23) but with a slightly different damping function:
fCHGdamp(RAB)=[1+a(RAB/R0,AB)-12]-1 | (5.26) |
This version is activated by setting DFT_D = EMPIRICAL_CHG, and the damping parameter a is controlled by the keyword DFT_D_A.
DFT_D_A
DFT_D_A
Controls the strength of dispersion corrections in the Chai–Head-Gordon DFT-D scheme, Eq. (5.26).
TYPE:
INTEGER
DEFAULT:
600
OPTIONS:
n
Corresponding to a=n/100.
RECOMMENDATION:
Use the default.
Note: 1. DFT-CHG is only defined for elements up to Xe. 2. The ωB97X-D and ωM05-D functionals automatically employ the DFT-CHG dispersion correction when used via METHOD = wB97X-D or wM05-D.
Grimme’s DFT-D3 method
450
J. Chem. Phys.
(2010),
132,
pp. 154104.
Link
constitutes an improvement
on his D2 approach, and is also available along with analytic first and second derivatives,
for any density functional that is available in Q-Chem. The D3 correction
includes a potential akin to that in D2 but including atomic C8 terms as well:
ED3,2-body=-atoms∑Aatoms∑B<A[s6(C6,ABR6AB)fdamp,6(RAB)+s8(C8,ABR8AB)fdamp,8(RAB)]. | (5.27) |
The total D3 dispersion correction consists of this plus a three-body term of the Axilrod-Teller-Muto (ATM) triple-dipole variety, so that the total D3 energy is EDFT-D3=EKS-DFT+ED3,2-body+EATM,3-body
Several versions of DFT-D3 are available as of Q-Chem 5.0, which differ in
the choice of the two damping functions. Grimme’s
formulation,
450
J. Chem. Phys.
(2010),
132,
pp. 154104.
Link
which is now known as the “zero-damping”
version [DFT-D3(0)], uses damping functions of the form
fD3(0)damp,n(RAB)=[1+6(RABsr,nR0,AB)-βn]-1 | (5.28) |
for n=6 or 8, β6=14, and β8=16.
450
J. Chem. Phys.
(2010),
132,
pp. 154104.
Link
,
786
J. Chem. Theory Comput.
(2013),
9,
pp. 263.
Link
The parameters R0,AB
come from atomic van der Waals radii, sr,6 is a functional-dependent
parameter, and sr,8=1. Typically s6 is set to unity and s8 is
optimized for the functional in question.
The more recent Becke–Johnson-damping version of DFT-D3,
452
J. Comput. Chem.
(2011),
32,
pp. 1456.
Link
DFT-D3(BJ), is designed to be finite (but non-zero) as RAB→0.
The damping functions used in DFT-D3(BJ) are
fD3(BJ)damp,n(RAB)=RnABRnAB+(α1R0,AB+α2)n | (5.29) |
where α1 and α2 are adjustable parameters fit for each
density functional. As in DFT-D3(0), s6 is generally fixed to unity and
s8 is optimized for each functional. DFT-D3(BJ) generally outperforms the
original DFT-D3(0) version.
452
J. Comput. Chem.
(2011),
32,
pp. 1456.
Link
The C6-only (CSO) approach of Schröder et al.
1136
J. Chem. Theory Comput.
(2015),
11,
pp. 3163.
Link
discards
the C8 term in Eq. (5.27) and uses a damping function with one
parameter, α1:
fD3(CSO)damp,6(RAB)=C6ABR6AB+(2.5Å)6(s6+α11+exp[RAB-(2.5Å)R0,AB]). | (5.30) |
The DFT-D3(BJ) approach was re-parameterized by Smith et al.
1190
J. Phys. Chem. Lett.
(2016),
7,
pp. 2197.
Link
to
yield the “modified” DFT-D3(BJ) approach, DFT-D3M(BJ), whose parameterization
relied heavily on non-equilibrium geometries. The same authors also introduces
a modification DFT-D3M(0) of the original zero-damping correction, which
introduces one additional parameter
(α1) as compared to DFT-D3(0):
fD3M(0)damp,n(RAB)=[1+6(RABsr,nR0,AB+α1R0,AB)-βn]-1. | (5.31) |
Finally, optimized power approach of Witte et al.
1378
J. Chem. Theory Comput.
(2017),
13,
pp. 2043.
Link
treats the
exponent, β6, as an optimizable parameter, given by
fD3(op)damp,n(RAB)=RβnABRβnAB+(α1R0,AB+α2)βn. | (5.32) |
Note that β8=β6+2.
To summarize this bewildering array of D3 damping functions:
DFT-D3(0) is requested by setting DFT_D = D3_ZERO. The model depends on four scaling parameters (s6, sr,6, s8, and sr,8), as defined in Eq. (5.28).
DFT-D3(BJ) is requested by setting DFT_D = D3_BJ. The model depends on four scaling parameters (s6, s8, α1, and α2), as defined in Eq. (5.29).
DFT-D3(CSO) is requested by setting DFT_D = D3_CSO. The model depends on two scaling parameters (s6 and α1), as defined in Eq. (5.30).
DFT-D3M(0) is requested by setting DFT_D = D3_ZEROM. The model depends on five scaling parameters (s6, s8, sr,6, sr,8, and α1), as defined in Eq. (5.31).
DFT-D3M(BJ) is requested by setting DFT_D = D3_BJM. The model depends on four scaling parameters (s6, s8, α1, and α2), as defined in Eq. (5.29).
DFT-D3(op) is requested by setting DFT_D = D3_OP. The model depends on four scaling parameters (s6, s8, α1, α2, and β6), as defined in Eq. (5.29).
The scaling parameters in these damping functions can be modified using the $rem variables described below. Alternatively, one may simply set DFT_D = D3, and a D3 dispersion correction will be selected automatically, if one is available for the selected functional.
Note: 1. DFT-D3(0) is defined for elements up to Pu (Z=94). 2. The B97-D3(0), ωB97X-D3, ωM06-D3 functionals automatically employ the DFT-D3(0) dispersion correction when invoked by setting METHOD equal to B97-D3, wB97X-D3, or wM06-D3. 3. The old way of invoking DFT-D3, namely through the use of EMPIRICAL_GRIMME3, is still supported, though its use is discouraged since D3_ZERO accomplishes the same thing but with additional precision for the relevant parameters. 4. When DFT_D = D3, parameters may not be overwritten, with the exception of DFT_D3_3BODY; this is intended as a user-friendly option. This is also the case when EMPIRICAL_GRIMME3 is employed for a functional parameterized in Q-Chem. When any of D3_ZERO, D3_BJ, etc. are chosen, Q-Chem will automatically populate the parameters with their default values, if they available for the desired functional, but these defaults can still be overwritten by the user.
DFT_D3_S6
DFT_D3_S6
The linear parameter s6 in eq. (5.27). Used in all forms of DFT-D3.
TYPE:
INTEGER
DEFAULT:
100000
OPTIONS:
n
Corresponding to s6=n/100000.
RECOMMENDATION:
NONE
DFT_D3_RS6
DFT_D3_S8
DFT_D3_S8
The linear parameter s8 in Eq. (5.27). Used in DFT-D3(0),
DFT-D3(BJ), DFT-D3M(0), DFT-D3M(BJ), and DFT-D3(op).
TYPE:
INTEGER
DEFAULT:
100000
OPTIONS:
n
Corresponding to s8=n/100000.
RECOMMENDATION:
NONE
DFT_D3_RS8
DFT_D3_A1
DFT_D3_A2
DFT_D3_POWER
DFT_D3_POWER
The nonlinear parameter β6 in Eq. (5.32). Used in
DFT-D3(op). Must be greater than or equal to 6 to avoid divergence.
TYPE:
INTEGER
DEFAULT:
600000
OPTIONS:
n
Corresponding to β6=n/100000.
RECOMMENDATION:
NONE
The three-body interaction term, E(3),
450
J. Chem. Phys.
(2010),
132,
pp. 154104.
Link
must be
explicitly turned on, if desired.
DFT_D3_3BODY
DFT_D3_3BODY
Controls whether the three-body interaction in Grimme’s DFT-D3 method should be applied
(see Eq. (14) in Ref.
450
J. Chem. Phys.
(2010),
132,
pp. 154104.
Link
).
TYPE:
LOGICAL
DEFAULT:
FALSE
OPTIONS:
FALSE (or 0)
Do not apply the three-body interaction term
TRUE
Apply the three-body interaction term
RECOMMENDATION:
NONE
More recently, Grimme published an extended D3 model,
D4.
162
J. Chem. Phys.
(2017),
147,
pp. 034112.
Link
,
163
J. Chem. Phys.
(2019),
150,
pp. 154122.
Link
,
164
Phys. Chem. Chem. Phys.
(2020),
22,
pp. 8499.
Link
The main feature
of D4 is that the coefficients are generated through Casimir-Polder integration
of the dynamic atomic polarizabilities α(iω) where electronic
density information is employed via atomic partial charges. Benchmark results
show that the proposed D4 model yields significantly lower error bars. The
DFT-D4 dispersion energy similar to D3 model is given by
ED4,2-body=-atoms∑Aatoms∑B<A[s6(C6,ABR6AB)fBJ damp,6(RAB)+s8(C8,ABR8AB)fBJ damp,8(RAB)]. | (5.33) |
The Becke-Johnson damping is utilized as default. The coordination number dependent CAB6 coefficients are obtained on-the-fly via Casimir-Polder integration
CAB6=NA,ref∑A,ref=1NB,ref∑B,ref=13π∫∞0𝑑ωαA,ref(iω,zA)×WA,refAαB,ref(iω,zB)WB,refB, | (5.34) |
where
αA,ref(iω,zA)=1m[αAmHn(iω)-n2αH2(iω)×ζ(zHA,ref,zH2)]ζ(zA,zA,ref) | (5.35) |
and
ζ(zA,zA,ref)=ba[1.47exp(zA/zA,ref)log10(zA,ref/zA)]. | (5.36) |
αAmHn denotes the reference polarizailities which represents the molecular polarizability of symmetric hydride systems AmHn. WA,ref/B,refA/B are weighting factors determining the contributions of all element specific reference systems NA,ref/B,ref. zHA describes the effective charge of hydrogen connectd to atom A in the reference system AmHn. The effective charge zA is computed self-consistently via Mulliken charge qA,
zA=ZA+qA. | (5.37) |
The coefficients a and b are parametrized to match cationic static polarzizabilies and TD-DFT derived molecular dispersion coefficients, respectively.
DFT_D4_S6
DFT_D4_S6
The linear parameter s6. Used in DFT-D4.
TYPE:
INTEGER
DEFAULT:
Optimized number for the specified functional
OPTIONS:
n
Corresponding to s6=n/100000000.
RECOMMENDATION:
NONE
DFT_D4_S8
DFT_D4_S8
The linear parameter s8. Used in DFT-D4.
TYPE:
INTEGER
DEFAULT:
Optimized number for the specified functional
OPTIONS:
n
Corresponding to s8=n/100000000.
RECOMMENDATION:
NONE
DFT_D4_S10
DFT_D4_S10
The linear parameter s10. Used in DFT-D4.
TYPE:
INTEGER
DEFAULT:
Optimized number for the specified functional
OPTIONS:
n
Corresponding to s10=n/100000000.
RECOMMENDATION:
NONE
DFT_D4_A1
DFT_D4_A1
The nonlinear parameter α1. Used in DFT-D4.
TYPE:
INTEGER
DEFAULT:
Optimized number for the specified functional
OPTIONS:
n
Corresponding to α1=n/100000000.
RECOMMENDATION:
NONE
DFT_D4_A2
DFT_D4_A2
The nonlinear parameter α2. Used in DFT-D4.
TYPE:
INTEGER
DEFAULT:
Optimized number for the specified functional
OPTIONS:
n
Corresponding to α2=n/100000000.
RECOMMENDATION:
NONE
DFT_D4_S9
DFT_D4_S9
The linear parameter s9. Used in DFT-D4.
TYPE:
INTEGER
DEFAULT:
Optimized number for the specified functional
OPTIONS:
n
Corresponding to s9=n/100000000.
RECOMMENDATION:
NONE
DFT_D4_WF
DFT_D4_WF
Weighting factor for Gaussian weighting.
TYPE:
INTEGER
DEFAULT:
600000000
OPTIONS:
n
Corresponding to wf=n/100000000.
RECOMMENDATION:
Use default
DFT_D4_GA
DFT_D4_GA
Charge scaling
TYPE:
INTEGER
DEFAULT:
300000000
OPTIONS:
n
Corresponding to ga=n/100000000.
RECOMMENDATION:
Use default
DFT_D4_GC
DFT_D4_GC
Charge scaling
TYPE:
INTEGER
DEFAULT:
200000000
OPTIONS:
n
Corresponding to gc=n/100000000.
RECOMMENDATION:
Use default
Example 5.11 Applications of B3LYP-D3(0) with custom parameters to a methane dimer.
$comment Geometry optimization, followed by single-point calculations using a larger basis set. $end $molecule 0 1 C 0.000000 -0.000323 1.755803 H -0.887097 0.510784 1.390695 H 0.887097 0.510784 1.390695 H 0.000000 -1.024959 1.393014 H 0.000000 0.001084 2.842908 C 0.000000 0.000323 -1.755803 H 0.000000 -0.001084 -2.842908 H -0.887097 -0.510784 -1.390695 H 0.887097 -0.510784 -1.390695 H 0.000000 1.024959 -1.393014 $end $rem JOBTYPE opt EXCHANGE B3LYP BASIS 6-31G* DFT_D D3_ZERO DFT_D3_S6 100000 DFT_D3_RS6 126100 DFT_D3_S8 170300 DFT_D3_3BODY FALSE $end @@@ $molecule read $end $rem JOBTYPE sp EXCHANGE B3LYP BASIS 6-311++G** DFT_D D3_ZERO DFT_D3_S6 100000 DFT_D3_RS6 126100 DFT_D3_S8 170300 DFT_D3_3BODY FALSE $end