X

Search Results

Searching....

5.7 DFT Methods for van der Waals Interactions

5.7.3 Empirical Dispersion Corrections: DFT-D

(February 4, 2022)

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.2; 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).400, 395

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 Grimme399 EMPIRICAL_CHG DFT-CHG dispersion correction from Chai and Head-Gordon183 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.392 D3_BJ DFT-D3(BJ) dispersion correction from Grimme et al.394 D3_CSO DFT-D3(CSO) dispersion correction from Schröder et al.1006 D3_ZEROM DFT-D3M(0) dispersion correction from Smith et al.1051 D3_BJM DFT-D3M(BJ) dispersion correction from Smith et al.1051 D3_OP DFT-D3(op) dispersion correction from Witte et al.1226 D3 Automatically select the “best” available D3 dispersion correction D4 DFT-D4 dispersion correction from Caldeweyher et al.146, 147, 148
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,399 in which the empirical dispersion potential has the aforementioned form, namely, pairwise atomic -C/R6 terms:

EdispD2=-s6AatomsB<Aatoms(C6,ABRAB6)fdampD2(RAB). (5.23)

This function is damped at short range, where RAB-6 diverges, via

fdampD2(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+EdispD2.

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.2 seems inconsistent and is not recommended. The global parameter s6 in Eq. (5.23) was optimized by Grimme for four different functionals,399 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,183 which uses the same form as Eq. (5.23) but with a slightly different damping function:

fdampCHG(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 method392 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=-AatomsB<Aatoms[s6(C6,ABRAB6)fdamp,6(RAB)+s8(C8,ABRAB8)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,392 which is now known as the “zero-damping” version [DFT-D3(0)], uses damping functions of the form

fdamp,nD3(0)(RAB)=[1+6(RABsr,nR0,AB)-βn]-1 (5.28)

for n=6 or 8, β6=12, and β8=14. 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,394 DFT-D3(BJ), is designed to be finite (but non-zero) as RAB0. The damping functions used in DFT-D3(BJ) are

fdamp,nD3(BJ)(RAB)=RABnRABn+(α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.394

The C6-only (CSO) approach of Schröder et al.1006 discards the C8 term in Eq. (5.27) and uses a damping function with one parameter, α1:

fdamp,6D3(CSO)(RAB)=CAB6RAB6+(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.1051 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):

fdamp,nD3M(0)(RAB)=[1+6(RABsr,nR0,AB+α1R0,AB)-βn]-1. (5.31)

Finally, optimized power approach of Witte et al.1226 treats the exponent, β6, as an optimizable parameter, given by

fdamp,nD3(op)(RAB)=RABβnRABβn+(α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_RS6
       The nonlinear parameter sr,6 in Eqs. (5.28) and Eq. (5.31). Used in DFT-D3(0) and DFT-D3M(0).
TYPE:
       INTEGER
DEFAULT:
       100000
OPTIONS:
       n Corresponding to sr,6=n/100000.
RECOMMENDATION:
       NONE

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_RS8
       The nonlinear parameter sr,8 in Eqs. (5.28) and Eq. (5.31). Used in DFT-D3(0) and DFT-D3M(0).
TYPE:
       INTEGER
DEFAULT:
       100000
OPTIONS:
       n Corresponding to sr,8=n/100000.
RECOMMENDATION:
       NONE

DFT_D3_A1

DFT_D3_A1
       The nonlinear parameter α1 in Eqs. (5.29), (5.30), (5.31), and (5.32). Used in DFT-D3(BJ), DFT-D3(CSO), DFT-D3M(0), DFT-D3M(BJ), and DFT-D3(op).
TYPE:
       INTEGER
DEFAULT:
       100000
OPTIONS:
       n Corresponding to α1=n/100000.
RECOMMENDATION:
       NONE

DFT_D3_A2

DFT_D3_A2
       The nonlinear parameter α2 in Eqs. (5.29) and (5.32). Used in DFT-D3(BJ), DFT-D3M(BJ), and DFT-D3(op).
TYPE:
       INTEGER
DEFAULT:
       100000
OPTIONS:
       n Corresponding to α2=n/100000.
RECOMMENDATION:
       NONE

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),392 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. 392).
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.146, 147, 148 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=-AatomsB<Aatoms[s6(C6,ABRAB6)fBJ damp,6(RAB)+s8(C8,ABRAB8)fBJ damp,8(RAB)]. (5.33)

The Becke-Johnson damping is utilized as default. The coordination number dependent C6AB coefficients are obtained on-the-fly via Casimir-Polder integration

C6AB=A,ref=1NA,refB,ref=1NB,ref3π0𝑑ωαA,ref(iω,zA)×WAA,refαB,ref(iω,zB)WBB,ref, (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/BA,ref/B,ref 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

View output