SAPT(DFT) is an extension of the SAPT theory. In SAPT(DFT), the monomers are described in terms of Kohn-Sham (KS) orbitals and orbital energies as well as of TD-DFT response functions. A complete description of the theory together with references to the historical developments of the method, as well as numerical examples, can be found in Ref. 14.
SAPT(DFT) consists of two steps. First is the so-called SAPT(KS), where electrostatics, first-order exchange, induction, exchange-induction, dispersion, and exchange-dispersion are obtained by using SAPT terms of the zeroth order in W with SCF orbitals and orbital energies replaced by their KS equivalents. For meaningful results, SAPT(KS) requires asymptotically corrected (AC) Kohn-Sham calculations. SAPT(KS) does not reproduce dispersion correctly. This problem was found to be due to the use of a formula asymptotically related to uncoupled dynamic polarizabilities. Instead, the dispersion (and induction) energies should be calculated from frequency-dependent density susceptibility (FDDS) functions, also referred to as propagators, obtained from TD-DFT, i.e., at the coupled Kohn-Sham (CKS) level. The SAPT(KS) and CKS steps form the complete SAPT(DFT). The total SAPT(DFT) interaction energy (up to second order in V ) can be defined as [14]
where the terms with the CKS label result from the coupled Kohn–Sham approach. The Ẽ_{exch-ind}^{(2)}(CKS) and Ẽ_{exch-disp}^{(2)}(CKS) are approximated by scaling their KS counterparts as described in Ref. 14. The method can be shown to be potentially exact for all major components of the interaction energy (asymptotically for exchange interactions) in the sense that these components would be exact if the DFT description of the monomers were exact.The nominal scaling of SAPT(KS) is O(N^{5}) and that of SAPT(DFT) is O(N^{6}), in both cases significantly better than the O(N^{7}) scaling of the regular SAPT. The scaling of SAPT(DFT) can be lowered to O(N^{5}) by using density fitting [12, 14, 15, 57, 62]. Furthermore, although the density-fitting transformation still scales as O(N^{5}), it has a greatly reduced prefactor.
SAPT(DFT) requires dalton
2.0 [37] or orca [53] for monomer DFT calculations (in the
latter case, only the density-fitting version has been interfaced, because
orca does not provide two-electron integrals needed
for the standard version. See Sec. 16.4.2
for the detailed description). The
users of SAPT(DFT) need to obtain a separate license and download the
package from
http://www.kjemi.uio.no/software/dalton/dalton.html.
Consult the dalton manual for
installation requirements. Some of the SAPT(DFT) functionality is available
via a patch dalton.diff distributed together
with sapt2012. First, compile and install the standard
dalton 2.0. After testing the installation, update
the path to the dalton directory in
the patchdalton script located in the main
sapt2012 directory and run the script. If the
patching is successful, compile the patched dalton. The recompilation should last significantly
shorter than the original compilation. If you have not compiled
sapt2012 with dalton support, set
the DALTON variable of the Compall script to the proper path to the dalton executable, make sure the variable
SAPTDFT=YES, and recompile sapt2012. If the compilation is error-free, the
SAPT(DFT) code is ready to use (both the regular and density-fitted
versions are created).
Before performing calculations for systems of interest to you, check the examples from examples/saptdft directory. The SAPT(DFT) run consists of the monomers calculation, transformation, the SAPT(KS) step, and the CKS step. Notice that some of the errors printed by Dalton are harmless in the SAPT(DFT) calculations. In particular,
is a normal stop after calculating the necessary one-electron integrals. Warning. Patched dalton should not be used for any ather purpose except for SCF/DFT calculations. The effect of the patch on other parts of the dalton code has not been tested.
The keywords discussed below are for both the regular and density-fitted versions of SAPT(DFT). Some additional keywords for the latter case are described in Sec. 16.4 . SAPT(DFT) requires SAPTKS=T keyword in the INPUTCOR namelist of the nameP.data file and proper keywords in the TRN namelist prepared like for the regular runs (see also Sec. 9.5 ). For CKS calculations that are a part of SAPT(DFT), a separate SAPTDFT namelist with keywords CKSDISP=T and CKSIND=T should also be present at the minimum. The complete set of SAPTDFT namelist keywords is:
The dalton patch enables some new options in dalton which are required for SAPT(DFT) calculations. All the new keywords belong to the *DFT INPUT submodule. The following new keywords are implemented:
For the DFT calculations, we recommend the standard Dalton grid or one of the denser grids. For the .THRESH parameter in *HF INPUT, we recommend values of about 10^{-7}.
Since SAPT(DFT) currently works with dalton interface only (and with orca in the density-fitted version), a few scripts are located in misc/daltutil to help in creating dalton input files. Script atmol2dalton converts atmol name.intinp input files into name.mol required by dalton, altergeo.awk updates the name.mol geometry, and with createinputs one can construct name.mol files using basis sets included in the dalton distribution. The usage of the scripts is explained in the sources. The name.dal files can be taken from the examples.
The SAPT(DFT) formalism presented above has been implemented up to second order in the operator V . For systems with a significant induction component (e.g., the water dimer), some higher-order terms constitute a large fraction of the total interaction energy. Those terms can be estimated from the supermolecular HF-SCF approach, with δE_{int,resp}^{HF} defined in Eq. (5 ). Since no such term can be computed from a pure DFT calculation, this term has to be extracted from a separate SAPT/HF-SCF run. For the calculation of δE_{int,resp}^{HF}, no expensive SAPT terms are needed and the scaling is O(N^{5}). Nevertheless, the time of this calculation would be similar to that of the SAPT(DFT) step since the costly transformation step must be done again. To calculate only δE_{int,resp}^{HF}, one should prepare inputs as for the regular SAPT runs (any interface can be used as long as the geometry and basis sets are the same as in the SAPT(DFT) run), set DELTASCF=T in the INPUTCOR namelist (cf. Sec. 10.2.3 ), and use scfcp option for the SAPT script. The δE_{int,resp}^{HF} value should then be added to the SAPT(DFT) interaction energy. An example of such a calculation is in the directory examples/saptdft/H2O-aTZ-MC.
Density-fitting (DF) method (also known as resolution of identity) has been implemented in SAPT(DFT) [15, 57]. The resulting approach has the scaling reduced from O(N^{6}) to O(N^{5}). Substantial savings are also achieved for O(N^{5}) terms. The code with density fitting is also more memory and disk efficient.
Running DF-SAPT(DFT) requires some modifications of the input files and the creation of a file with an auxiliary basis set. The latter file should be named name.aux and contains coefficients of the auxiliary basis function for each atom. The format is the following (see also examples/dfsaptdft):
where charge is the charge of an atom and the coefficients are in the turbomole format. We recommend basis sets from Ref. 61. These basis sets are available in turbomole format at ftp://ftp.chemie.uni-karlsruhe.de/pub/cbasen/. The number and ordering of the atoms must be identical to the main basis set used. The quality of the auxiliary basis set is critical for the accuracy of the DF approach and using auxiliary basis set optimized for a different main basis would result in poor accuracy.
In the TRN namelist of the nameP.data, one should specify T2EL=F to suppress the standard two-electron transformation. There is also an additional DF namelist. The only parametr used in this memlist is MEMTRAN=x keyword, where x is equal to either the number of requested memory words or zero. In the latter case, the memory is allocated automatically to the lowest reasonable level. Finally, for running DF-SAPT(DFT), one uses SAPTdf and RunlotDALTONdf scripts instead of SAPT and RunlotDALTON, respectively. Example inputs and scripts for running DF-SAPT(DFT) (analogous to non-DF SAPT(DFT) counterparts in examples/saptdft) are located in examples/dfsaptdft, along with the pertinent outputs. Notice small DF errors when comparing the results of DF and non-DF calculations.
For larger auxiliary basis sets, where the linear dependencies start to emerge, the numerical errors start to amplify and the accuracy of the electrostatic energy is diminished. The problem has been described in Ref. 57. An effective solution is to perform the inversion of the auxiliary J matrix with quadruple precision (QP). For compilers supporting QP, one should use this option to improve the accuracy. The QP code has been tested with ifort and ibm64 platforms. To compile the QP code, one should include QUADINV=YES in the Compall script. If the compiler does not support quadruple precision, the compilation will fail. Since QP calculations are more time consuming than double precision ones, the QP inversion is not turned on by default. To use it, one needs to also specify quad option for SAPTdf or RunlotDALTONdf scripts. With this option, after the SAPT(DFT) calculation is complete, an extra calculation of the electrostatic energy is performed. This result replaces the standard double precision results in the summary table. By comparing the double and quadruple precision values, one can estimate the severity of the quasi-linear dependence. Additionally, slightly better accuracy of density-fitting of the electrostatic energy is obtained by using basis sets of Ref. 65 rather than those of Ref. [61]. One can include such basis set in name-elst.aux and it will be used instead of name.aux in the quadruple precision calculation of the electrostatic energy. If the file is missing, the standard name.aux will be used.
DF-SAPT(DFT) is also partially parallelized. All calculations performed by the (patched) dalton program can be run in parallel. Also the calculations of the CKS kernel integrals are parallelized. The rest of the programs is not parallelized but can be used with parallel dalton.
DF-SAPT(DFT) is limited to up to g functions in the auxiliary basis set. This restriction is due to the use of gamess integrals. In addition, when using the orca interface, the main basis set is currently limited to up to f functions (in practice, both these restrictions are equivalent, because the auxiliary basis set should contain functions with the highest angular momentum L at least by one larger than the main basis set).
The density-fitting code contains some atomic integrals code from gamess. The gamess authors requested that the users of density-fitting SAPT(DFT) should cite gamess [22] along with SAPT.
The dalton input files required in a DF-SAPT(DFT) run are almost identical with those in regular SAPT(DFT). However, in name.dal files, CKS keyword should be replaced by CKSAUX. This keyword works only with the LDA kernel, therefore, GGAKER keyword cannot be used. It it recommended to include .NOTWO keyword in **INTEGRALS of nameA.dal and nameB.dal when monomer basis sets are used.
Orca is a quantum chemistry package developed by Neese and coworkers [53]. Its DFT module has very efficient density-fitting techniques and we recommend this front-end for DF-SAPT(DFT) calculations, especially in the case of monomers larger than about a dozen of atoms. The Linux, Windows, and Mac-OS X binaries are available at http://www.thch.uni-bonn.de/tc/orca, free of charge for academic users (a registration and acquiring a license is necessary). The orca source files are not needed by DF-SAPT(DFT). Version 2.9.0 or later is required, because earlier versions contain errors in the GRAC asymptotic correction code and lead to incorrect SAPT results. The installation of orca can be done before or after the installation of SAPT, as both are completely independent, but in the latter case one has to edit the file vars.cfg so that the variable ORCADIR is set to the directory with the orca executables.
To use DF-SAPT(DFT) with orca, ISITORCA=T must be set in the TRN namelist of the nameP.data file. The required monomer and/or dimer orca input files follow the same naming pattern (ending with .data) as in the case of gaussian (see Sec. 10.1 for the detailed description). Each such file, besides containing regular orca input in the first section, must end with a special second section starting with a line containing the string META as its first four characters. This section always defines the molecular geometry and basis set used in the DFT calculation (even if a basis set is already defined in the first section, it gets overridden here). Additionally, if one of the density-fitting DFT algorithms implemented in orca is invoked, the META section can be used to override the auxiliary (density fitting) basis set defined before or provided by orca as a default. The example shown below can be used to performed the DFT calculation of the first monomer within the water dimer, using the full dimer-centered basis set (file nameA.data).
The first line declares the use of the PBE functional, def2-SVP as the main basis set and def2-SVP/J as the Coulomb-fitting auxiliary basis set, as well as tight convergence criteria for the SCF procedure and the use of bohr units in the geometry section (see the orca manual for details). The gracLB true setting is essential for any SAPT(DFT) calculation as it turns on the GRAC asymptotic correction of the exchange-correlation functional (the other type of asymptotic correction, Fermi-Amaldi, is curently not implemented in orca). The correction requires the ionization potential of the molecule, which must be given in electron-volt units (as opposed to the hartree units used in dalton). The META section starts with a line containing three integers: number of atoms (including ghosts), total charge and multiplicity. For each atom, there is a line containing: atomic number, x, y, and z Cartesian coordinates, and atomic charge (zero for ghost atoms). The main basis set for each element is defined in a block (separated by exactly one empty line from previous blocks) starting with a line containing the string BAS as the first three characters and the atomic number. The next line of the BAS block gives the total number of contracted orbitals. For each such orbital, there is a sub-block starting with two integers (angular momentum and number of primitive functions) and then containing, in separates lines, pairs of exponential Gaussian parameters and contraction coefficients. The order of elements (BAS blocks) is arbitrary but all elements contained in the geometry specification must have an explicitly defined basis set (in particular, it can be specified as empty by setting the number of contracted orbitals to zero), even if this definition is redundant with respect to the basis set keyword in the general orca input. This requirement stems from the fact that the META section information is used not only by orca but also by a separate one-electron integral module, based on routines extracted from the gamess(us) code. Currently, only Gaussian basis functions with angular momenta up to f can be used.
The def2-SVP/J auxiliary basis set requested in the orca input could be replaced by another basis set, by using in the META section additional blocks, with the keyword BAS replaced by AUX. The order of these blocks is also arbitrary but no empty AUX blocks are allowed (this is an orca restriction). Note also that the RI-J Coulomb fitting (requested by the RI keyword in our example) can only be used with pure density functionals. For hybrid functionals, orca has two other methods known as RIJCOSX and RI-JK, the latter requiring “/JK” (rather than “/J”) auxiliary basis sets. All these DFT fitting basis sets must not be confused with auxiliary basis sets required by DF-SAPT(DFT) and defined in the file name.aux, which are designed to fit the correlation effects and must be provided by the user, as in the case of dalton (such basis sets are actually used in orca at the MP2 level and denoted by “/C”).
The orca route in DF-SAPT(DFT) requires a preparation of one input file not occuring in the dalton route, called kernel.data. This file contains just two lines of the form
where n_rad and n_aug are the numbers of radial and angular grid points around each atom in the quadrature of the CKS kernel integral. In general, much looser grids can be used here compared to the main DFT grids without any significant loss of accuracy. Typically, the values of n_rad = 50 and n_ang = 50 lead to completely negligible errors, while values of n_rad = 30 and n_ang = 26 introduce errors in the dispersion energy of just a few tenths of one percent. Parameter mem specifies the maximum amount of memory (in 8-byte words) to be used in the kernel generation. For larger systems, limiting memory use requires slower multipass calculations. The parameter ksi is the fraction of exact exchange of the employed functional (for instance, 0.25 for PBE0). As the last difference compared to the dalton route, the script dfSAPT.orca rather than dfSAPT is used.