6.8 Dispersion Correction for DFT Calculations

Grimme DFT+Dispersion, DFT-D3, DFT-D4

Based on an idea that has earlier been proposed for Hartree-Fock calculations[192, 193], a general empirical dispersion correction has been proposed by Stefan Grimme for density functional calculations [194]. A modified version of the approach with extension to more elements and more functionals has been published in ref. [195]. A more recent implementation [196] is less empirical, i.e. the most important parameters are computed by first principles, and it provides a consistent description across the whole periodic system. This version, named DFT-D3 with Becke-Johnson damping, has been the default for dispersion correction ever since.

In 2017 a new model, termed D4, was published [197, 198] and released 2018. For details see below.

  • The first version (DFT-D1) can be invoked by the keyword $olddisp in the control file.

  • The second version (DFT-D2) is used if the keyword $disp is found.

  • For the usage of DFT-D3 just add keyword $disp3 to the control file. For DFT-D3 with Becke-Johnson damping, please add $disp3 -bj.

  • The latest version, DFT-D4, can be invoked by the keyword $disp4.

Only one of the four keywords is expected to be present.

DFT-D3

If DFT-D3 is used, the total energy is given by \[\begin{equation} E_{DFT-D3}=E_{KS-DFT}-E_{disp} \end{equation}\](6.18) where \(E_{KS-DFT}\) is the usual self-consistent Kohn-Sham energy as obtained from the chosen functional and \(E_{disp}\) is a dispersion correction given by the sum of two- and three-body energies \[\begin{equation} E_{disp}=E^{(2)}+E^{(3)}, \end{equation}\](6.19) with the dominating two-body term \[\begin{equation} E^{(2)}=\sum_{AB}\sum_{n=6,8,10,...} s_n\frac{C_n^{AB}}{r_{AB}^n} f_{d,n}(r_{AB}). \end{equation}\](6.20)

The first sum runs over all atom pair, \(C_n^{AB}\) denotes the \(n\)th-order dispersion coefficient for atom pair \(AB\), \(r_{AB}\) is their interatomic distance, and \(f_{d,n}\) is a damping function.

Becke-Johnson (BJ) damping can be invoked by adding the option bj or -bj to the $disp3 keyword: $disp bj If you use this damping option please also cite [199].

The three-body term can be switched on by adding abc to the $disp3 input line, i.e. to use it in combination with Becke-Johnson damping just add $disp3 bj abc

It is also possible not to use the functional name given in the control file but to tell the DFT-D3 routines to use the parameters which have been fitted to a specific functional. Just as in the original DFT-D3 routines, this can be selected by adding the func option, for example $disp3 bj func pbe0. It is recommended to use this option as the last one in the $disp3 input line.

Please have look at DFT-D3 homepage, Grimme group Bonn for more detailed information.

DFT-D4

DFT-D4 computes molecular dipole-dipole and dipole-quadrupole dispersion coefficients based on dynamic dipole polarizabilities calculated with TD-PBE38/daug-def2-QZVP used as atomic reference polarizabilities for elements up to radon (Z=86). In DFT-D4, all atomic reference polarizabilities are scaled according to Mulliken partial charges which are obtained by semi-empirical quantum mechanical tight binding, within the 2018 developed GFN2-xTB framework[71]. Alternatively a fallback option is provided in case GFN2-xTB fails to converge the electronic structure. Either classical, non-interative EN and CN dependent Gasteiger charges (option ’gasteiger’) or (externally provided) Hirshfeld-charges optained at least at PBE0/def2-TZVP level of theory can be used (option ’hirshfeld’).

For the scaling of the atomic reference polarizabilities, a special charge-function was designed containing one global parameter and the chemical hardnesses as elementspecific parameter to get a smooth scaling behavior. The charge-scaled dipole polarizabilities are Gaussian interpolated according to an empirically generated fractional coordination number to obtain the fully weighted dynamic polarizabilities \(\alpha(i\omega)\) for all atoms which are then numerically integrated via the Casimir-Polder scheme to obtain charge- and geometry-dependent dipole-dipole dispersion coefficients. The two-body energy expression has the usual sum over pair interactions form for dipole-dipole and dipole-quadrupole interactions.

Based on this expression, a self-consistent dispersion potential has been developed and implemented into the GFN2-xTB Hamiltonian to circumvent a costly coupled-perturbed SCF procedure when calculating analytical gradients. Furthermore, dynamic polarizabilities \(\alpha(i\omega)\) are used within an RPA-like expression to capture many-body interactions beyond the two-body terms.

The DFT-D4 energy expression differs for single point energies and gradient calculations, for reasons of computational efficiency. The gradient of the RPA-like dispersion energy scales as \(O(N^4)\) and cannot be screened, therefore the gradient is approximated by the ATM three-body dispersion expression as used in DFT-D3, which is computational less demanding and yields approximately similar gradients.

Parameters for the rational damping scheme introduced by Becke and Johnson are included for a large number of density functionals (including all standard functionals implemented in Turbomole) in the DFT-D4 code. The rational damping is the only scheme available in DFT-D4, also many-body dispersion terms are included by default.

The DFT-D4 keyword $disp4 can have, similar to $disp3, a couple of options. To see the list of options, just add -help to the keyword and inspect the output file.

As mentioned above, the gradient of the many body dispersion energy is done by using the Axilrod-Teller-Muto (ATM) three-body term, which is only an approximation compared to the term used in the energy evaluation. Hence, for larger molecules, the DFT-D4 energy and the DFT-D4 gradient do not match exactly. Very tight convergence criteria for both energy and gradient norm can not be expected in all cases.