SAPT: Symmetry-Adapted Perturbation Theory¶
Code author: Edward G. Hohenstein, Rob M. Parrish, Jérôme F. Gonthier, and Daniel. G. A. Smith
Section author: Edward G. Hohenstein and Jérôme F. Gonthier
Module: Keywords, PSI Variables, LIBSAPT_SOLVER
Warning
In rare cases with systems having a high degree of symmetry, PSI4 gives (very obviously) wrong answers for SAPT computations when the specification is in Z-matrix format. Use a Cartesian representation to avoid this problem.
Caution
In early versions (notably PSI4 alpha circa 2011 and before), frozen core was implemented incompletely and for only selected terms. Comparisons with papers published using early PSI4 SAPT code may show discrepancies of 0.01-0.10 kcal/mol in individual terms, particularly \(E_{exch}^{(11)}\) and \(E_{exch}^{(12)}\).
Caution
January 28th 2016, the default for all NAT_ORBS options was changed to true. Hence the code now by default uses natural orbital truncation to speed up the evaluation of energy terms wherever possible, according to literature recommendations. In early July 2016, some total SAPT energy psivars were renamed.
Caution
February 7, 2020, a missing term in \(E^{(30)}_{ind}\) was added, causing possible discrepancies with prior versions of the code on the order of 0.01 kcal/mol. See https://github.com/psi4/psi4/issues/1677
Caution
August 2021, the number of frozen core orbitals used in the dMP2 computations is now standardized. Specifically, we now rigorously enforce that the number of core orbitals frozen in dimer computations is equal to the sum of frozen orbitals of each monomer. Prior to this, a discrepency between these values was possible when one of the monomers was (exclusively) a charged alkali metal.
Caution
May 2022 c. v1.6, the default for DF_BASIS_ELST
changed from the value of DF_BASIS_SAPT (which itself
defaults to the RI of the orbital basis) to the JKFIT of the orbital
basis. This affects SAPT0 and sSAPT0 computed with the SAPT
module (the default code for energy("sapt0") that
can also compute higher-order SAPT). Electostatics, exchange,
and induction terms for SAPT0 and sSAPT0 accessed through
energy("sapt0") or energy("ssapt0") change; the dispersion
term does not change. The SAPT0 and sSAPT0 terms accessed as
subsidiary calculations of higher-order SAPT do not change; that is,
the SAPT module breaks the consistency of its SAPT0
results. The reasoning and reward behind this change is that the JKFIT
basis better describes the physics (see fitting changes ) and the
default SAPT0 results from the SAPT module are now
consistent with those from the FISAPT module and
the sapt(dft) module. See sapt-compare for an example.
To reproduce former behavior, set DF_BASIS_ELST to the
orbital basis set’s RI auxiliary basis.
Symmetry-adapted perturbation theory (SAPT) provides a means of directly computing the noncovalent interaction between two molecules, that is, the interaction energy is determined without computing the total energy of the monomers or dimer. In addition, SAPT provides a decomposition of the interaction energy into physically meaningful components: i.e., electrostatic, exchange, induction, and dispersion terms. In SAPT, the Hamiltonian of the dimer is partitioned into contributions from each monomer and the interaction.
Here, the Hamiltonian is written as a sum of the usual monomer Fock operators, \(F\), the fluctuation potential of each monomer, \(W\), and the interaction potential, \(V\). The monomer Fock operators, \(F_A+F_B\), are treated as the zeroth-order Hamiltonian and the interaction energy is evaluated through a perturbative expansion of \(V\), \(W_A\), and \(W_B\). Through first-order in \(V\), electrostatic and exchange interactions are included; induction and dispersion first appear at second-order in \(V\). For a complete description of SAPT, the reader is referred to the excellent review by Jeziorski, Moszynski, and Szalewicz [Jeziorski:1994:1887].
Several truncations of the closed-shell SAPT expansion are available in the SAPT module of PSI4. The simplest truncation of SAPT is denoted SAPT0 and defined in Eq. (1).
In this notation, \(E^{(vw)}\) defines the order in \(V\) and in \(W_A+W_B\); the subscript, \(resp\), indicates that orbital relaxation effects are included.
For SAPT(DFT), the SAPT expansion is truncated at the same level of SAPT0, but the zeroth-order Hamiltonian is chosen to be \(K_A+K_B\), the monomer Kohn-Sham operators, instead of the Fock operators. The dispersion term needs to be computed with orbital relaxation for the result to be reasonable, and this is possible by computing dispersion energy through coupled frequency-dependent density susceptibility (FDDS). The exchange-dispersion term is estimated by scaling the uncoupled exchange-dispersion energy.
The \(\delta_{HF}^{(2)}\) and \(\delta_{HF}^{(3)}\) terms take into account higher-order induction effects and are included in the definition of SAPT terms. They are computed from the Hartree–Fock supermolecular interaction energy \(E_{int}^{HF}\) and are only available in dimer-centered basis SAPT computations, which is the default (see below for monomer-centered basis computations). They are defined by:
Additionally, high-order coupling between induction and dispersion can be extracted from the supermolecular MP2 interaction energy:
where \(E_{int}^{MP2, corr}\) is the correlation part of the supermolecular MP2 interaction energy. \(\delta_{MP2}^{(2)}\) and \(\delta_{MP2}^{(3)}\) also improve the description of electrostatically dominated complexes. \(\delta_{MP2}^{(2)}\) can be applied to SAPT2+ or SAPT2+(3) energies whereas \(\delta_{MP2}^{(3)}\) should be applied to SAPT2+3 energies.
A thorough analysis of the performance of these truncations of closed-shell SAPT can be found in a review by Hohenstein and Sherrill [Hohenstein:2012:WIREs], and a systematic study of the accuracy of these truncations (with and without an improved CCD treatment of dispersion) using different basis sets is reported in [Parker:2014:094106].
The closed-shell SAPT module relies entirely on the density-fitting approximation of the two-electron integrals. The factorization of the SAPT energy expressions, as implemented in PSI4, assumes the use of density-fitted two-electron integrals, therefore, the closed-shell SAPT module cannot be run with exact integrals. In practice, we have found that the density-fitting approximation introduces negligible errors into the SAPT energy (often less than 0.01 kcal/mol for small dimers) and greatly improves efficiency.
The latest addition to the SAPT code is the SAPT0 method for open-shell
monomers [Gonthier:2016:134106]. This code is available for both exact and density fitted integrals,
except for the dispersion terms which implementation relies on a density fitting
factorization. Both UHF and ROHF REFERENCE can be used, but coupled
induction computations are currently not supported with ROHF. This means that orbital
relaxation is not included for ROHF and the uncoupled induction term is computed instead.
If both monomers are open-shell, their coupling is assumed to be high spin, i.e.
two doublets would interact to form a triplet.
The S2 approximation and scaling¶
All exchange terms in SAPT arise from the antisymmetrization of the wavefunctions of monomers A and B. Taking into account exchange of all possible electron pairs between the two monomers yields to complicated formulae. For this reason, exchange terms are often evaluated in the \(S^{2}\) approximation, that can be interpreted as the exchange of a single electron pair between monomers.
The \(S^{2}\) approximation is usually pretty good, but may break down for short intermolecular distance, particularly in high-order terms. To compensate for these deviations, in 2014 Parker et al. [Parker:2014:094106] recommended to scale all \(S^{2}\) approximated exchange terms by the ratio:
and recommended using the ratio with exponent \(\alpha = 1\). To obtain SAPT energies with this scaling,
simply set the keyword exch_scale_alpha true. Alternatively, another value for \(\alpha\)
can be specified by setting EXCH_SCALE_ALPHA to a value. For example,
set exch_scale_alpha 1.0
will set \(\alpha = 1.0\) and scale exchange energies with \(p_{EX}(1.0)\).
However, as pointed out by Schaffer and Jansen [Schaffer:2013:2570] in the context of DFT-based SAPT, the ratios \(E_{exch}^{(1)}(S^{2})/E_{exch}^{(1)}\), \(E_{\rm exch-ind}^{(2)}(S^2)/E_{\rm exch-ind}^{(2)}\), and \(E_{\rm exch-disp}^{(2)}(S^2)/E_{\rm exch-disp}^{(2)}\) are not very similar to each other. Hence, in 2016 universal scaling of all \(S^{2}\) approximated terms was turned off by default.
Finally, in 2014 Parker et al. [Parker:2014:094106] empirically discovered that SAPT0 energies for van der Waals dimers with close contacts benefit from a slightly modified recipe that involves an empirically adjusted exponent \(\alpha = 3.0\). To distinguish it from its unscaled counterpart, this energy is denoted sSAPT0 (see [Parker:2014:094106]).
where \(\delta_{HF}^{(2)}\) is computed without any scaling. Please note that
sSAPT0 is thus not the same as requesting exch_scale_alpha 3.0, and that the
scaling is automatically performed by requesting energy('ssapt0').
A First Example¶
The following is the simplest possible input that will perform all available SAPT computations (normally, you would pick one of these methods, not all of them).
molecule water_dimer {
0 1
O -1.551007 -0.114520 0.000000
H -1.934259 0.762503 0.000000
H -0.599677 0.040712 0.000000
--
0 1
O 1.350625 0.111469 0.000000
H 1.680398 -0.373741 -0.758561
H 1.680398 -0.373741 0.758561
units angstrom
no_reorient
symmetry c1
}
set basis aug-cc-pvdz
energy('sapt0')
energy('sapt2')
energy('sapt2+')
energy('sapt2+(3)')
energy('sapt2+3')
The SAPT module uses the standard PSI4 partitioning of the dimer
into monomers. SAPT does not use spatial symmetry and needs the geometry
of the system to remain fixed throughout monomer and dimer calculations.
These requirements are imposed whenever a SAPT calculation is requested
but can also be set explicitly with the no_reorient and symmetry
c1 molecule keywords, as in the example above. As a reminder, only
SAPT0 can handle the interaction of both closed- and open-shell monomers.
Higher-order SAPT and SAPT(DFT) are currently available only
for computation of interactions between
closed-shell singlets. The SAPT codes in PSI4 have been written
to utilize density fitting, which is much faster than using conventional
4-index electron repulsion integrals. This happens automatically and
does not require any additional keywords to be specified (both the
SAPT computations and the underlying Hartree-Fock computations will
utilize density fitting).
For SAPT(DFT), the user will need to manually specify the GRAC shift of both monomers, defined by the difference of ionization potential (IP) and Kohn-Sham HOMO energy. The ionization potential data for many common molecules is available in NIST Chemistry Webbook. Alternatively, one can estimate the ionization potential of molecule by computing the energy difference between the molecule as given, and the molecule after one electron has been removed (e.g., the energy difference between a neutral molecule and its cation).
The values of GRAC shifts should be given in Hartree. For the example above, the GRAC shift value of both molecules are calculated to be 0.1307 (using experimental IP from CCCBDB), and the user would specify them using the following keywords:
set globals {
sapt_dft_grac_shift_a 0.1307
sapt_dft_grac_shift_b 0.1307
}
A complete, minimal example of a SAPT(DFT) computation is given below.
molecule {
0 1
O -1.551007 -0.114520 0.000000
H -1.934259 0.762503 0.000000
H -0.599677 0.040712 0.000000
--
0 1
O 1.350625 0.111469 0.000000
H 1.680398 -0.373741 -0.758561
H 1.680398 -0.373741 0.758561
units angstrom
}
set {
basis aug-cc-pVDZ
sapt_dft_grac_shift_a 0.1307
sapt_dft_grac_shift_b 0.1307
}
energy('sapt(dft)')
As already mentioned above, the SAPT0 module for open-shell cases can also use exact integrals for all terms except for dispersion. In practice, density fitting is considerably faster and introduces negligible errors, thus it is the preferred method for open-shell cases as well. Below, you can find a minimum example of open-shell SAPT0 computation.
molecule {
0 1
O 0.000000 0.000000 6.000000
H 0.000000 1.431500 4.890600
H 0.000000 -1.431500 4.890600
--
0 2
O 0.000000 0.000000 0.000000
O 0.000000 2.503900 0.000000
H 0.000000 -0.424700 -1.839500
units bohr
symmetry c1
no_reorient
no_com
}
set {
reference uhf
scf_type df
basis aug-cc-pVDZ
}
energy('sapt0')
REFERENCE needs to be UHF or ROHF for the open-shell computation to proceed.
Advanced example¶
Open-shell computations can be difficult to converge in certain cases, thus you may want to have more control over the SCF procedure. You have the option of doing the driver job in the input file, by performing the dimer and monomer computations yourself. In the example below, we do a stability analysis for the open-shell monomer only
molecule {
0 2
O 0.000000 0.000000 0.000000
O 0.000000 2.503900 0.000000
H 0.000000 -0.424700 -1.839500
--
0 1
O 0.000000 0.000000 6.000000
H 0.000000 1.431500 4.890600
H 0.000000 -1.431500 4.890600
units bohr
symmetry c1
no_reorient
no_com
}
dimer = psi4.get_active_molecule()
set {
reference uhf
scf_type df
basis cc-pVDZ
df_basis_sapt cc-pVDZ-ri
guess sad
}
dimer = psi4.get_active_molecule()
set df_ints_io save
psi4.IO.set_default_namespace('dimer')
Edim, wfn_dimer = energy('scf',molecule=dimer,return_wfn=True)
set df_ints_io load
monomerA = dimer.extract_subsets(1,2)
psi4.IO.change_file_namespace(97, 'dimer', 'monomerA')
psi4.IO.set_default_namespace('monomerA')
set {
stability_analysis follow
}
EmonA, wfn_monA = energy('scf',molecule=monomerA,return_wfn=True)
monomerB = dimer.extract_subsets(2,1)
psi4.IO.change_file_namespace(97, 'monomerA', 'monomerB')
psi4.IO.set_default_namespace('monomerB')
set {
stability_analysis none
}
EmonB, wfn_monB = energy('scf',molecule=monomerB,return_wfn=True)
psi4.IO.change_file_namespace(97, 'monomerB', 'dimer')
psi4.IO.set_default_namespace('dimer')
aux_basis = psi4.core.BasisSet.build(wfn_dimer.molecule(), "DF_BASIS_SAPT",
psi4.core.get_global_option("DF_BASIS_SAPT"),
"RIFIT", psi4.core.get_global_option("BASIS"))
wfn_dimer.set_basisset("DF_BASIS_SAPT", aux_basis)
wfn_dimer.set_basisset("DF_BASIS_ELST", aux_basis)
psi4.sapt(wfn_dimer,wfn_monA,wfn_monB)
In this way, any of the SCF options can be tweaked for individual fragments.
For optimal efficiency, the example uses set df_ints_io save to keep file 97,
which contains the three-index integrals for density fitting. set df_ints_io load
then instructs the program to read these integrals from disk instead of recomputing
them. For each SCF computation, we use psi4.IO.set_default_namespace to uniquely
name scratch files. In the following SCF step, only file 97 is renamed using
psi4.IO.change_file_namespace so that integrals can be read from it.
For more information on stability analysis, see the stability
documentation.
SAPT0¶
Generally speaking, SAPT0 should be applied to large systems or large data sets. The performance of closed-shell SAPT0 relies entirely on error cancellation, which seems to be optimal with a truncated aug-cc-pVDZ basis, namely, jun-cc-pVDZ (which we have referred to in previous work as aug-cc-pVDZ’). We do not recommend using SAPT0 with large basis sets like aug-cc-pVTZ. A systematic study of the accuracy of closed-shell SAPT0 and other SAPT truncations, using different basis sets, is reported in [Parker:2014:094106]. In particular, an empirical recipe for scaled SAPT0 can yield improved performance and has been included in the output file as the sSAPT0 interaction energy. sSAPT0 is a free by-product and is automatically computed when SAPT0 is requested (see above for more details). The SAPT module has been used to perform SAPT0 computations with over 200 atoms and 2800 basis functions; this code should be scalable to 4000 basis functions. Publications resulting from the use of the SAPT0 code should cite the following publications: [Hohenstein:2010:184111] and [Hohenstein:2011:174107]. If the open-shell SAPT0 code is used, [Gonthier:2016:134106] should be additionally cited.
Basic SAPT0 Keywords¶
SAPT_LEVEL¶
The level of theory for SAPT
Type: string
Possible Values: SAPT0, SAPT2, SAPT2+, SAPT2+3
Default: SAPT0
BASIS¶
Primary basis set, describes the monomer molecular orbitals
Type: string
Possible Values: basis string
Default: No Default
DF_BASIS_SAPT¶
Auxiliary basis set for SAPT density fitting computations. Defaults to a RI basis.
Type: string
Possible Values: basis string
Default: No Default
DF_BASIS_ELST¶
Auxiliary basis set for SAPT Elst10 and Exch10 density fitting computations, may be important if heavier elements are involved. Defaults to a JKFIT basis. Previous to v1.6, defaulted to DF_BASIS_SAPT See fitting notes .
Type: string
Possible Values: basis string
Default: No Default
FREEZE_CORE¶
The scope of core orbitals to freeze in evaluation of SAPT \(E_{disp}^{(20)}\) and \(E_{exch-disp}^{(20)}\) terms. Recommended true for all SAPT computations
Type: string
Possible Values: FALSE, TRUE
Default: FALSE
D_CONVERGENCE¶
Convergence criterion for residual of the CPHF coefficients in the SAPT \(E_{ind,resp}^{(20)}\) term.
Type: conv double
Default: 1e-8
E_CONVERGENCE¶
Convergence criterion for energy (change) in the SAPT \(E_{ind,resp}^{(20)}\) term during solution of the CPHF equations.
Type: conv double
Default: 1e-10
MAXITER¶
Maximum number of CPHF iterations
Type: integer
Default: 50
PRINT¶
The amount of information to print to the output file for the sapt module. For 0, only the header and final results are printed. For 1, (recommended for large calculations) some intermediate quantities are also printed.
Type: integer
Default: 1
Advanced SAPT0 Keywords¶
AIO_CPHF¶
Do use asynchronous disk I/O in the solution of the CPHF equations? Use may speed up the computation slightly at the cost of spawning an additional thread.
Type: boolean
Default: false
AIO_DF_INTS¶
Do use asynchronous disk I/O in the formation of the DF integrals? Use may speed up the computation slightly at the cost of spawning an additional thread.
Type: boolean
Default: false
COUPLED_INDUCTION¶
Solve the CPHF equations to compute coupled induction and exchange-induction. These are not available for ROHF, and the option is automatically false in this case. In all other cases, coupled induction is strongly recommended. Only turn it off if the induction energy is not going to be used.
Type: boolean
Default: true
EXCH_SCALE_ALPHA¶
Whether or not to perform exchange scaling for SAPT exchange components. Default is false, i.e. no scaling. If set to true, performs scaling with \(Exch10 / Exch10(S^2)\). If set to a value \(\alpha\), performs scaling with \((Exch10 / Exch10(S^2))^{\alpha}\).
Type: string
Default: FALSE
INTS_TOLERANCE¶
Schwarz screening threshold. Minimum absolute value below which all three-index DF integrals and those contributing to four-index integrals are neglected. The default is conservative, but there isn’t much to be gained from loosening it, especially for higher-order SAPT.
Type: conv double
Default: 1.0e-12
DENOMINATOR_DELTA¶
Maximum error allowed (Max error norm in Delta tensor) in the approximate energy denominators employed for most of the \(E_{disp}^{(20)}\) and \(E_{exch-disp}^{(20)}\) evaluation.
Type: double
Default: 1.0e-6
DENOMINATOR_ALGORITHM¶
Denominator algorithm for PT methods. Laplace transformations are slightly more efficient.
Type: string
Possible Values: LAPLACE, CHOLESKY
Default: LAPLACE
DEBUG¶
The amount of information to print to the output file
Type: integer
Default: 0
Specific open-shell SAPT0 keywords¶
SAPT_MEM_SAFETY¶
Memory safety
Type: double
Default: 0.9
COUPLED_INDUCTION¶
Solve the CPHF equations to compute coupled induction and exchange-induction. These are not available for ROHF, and the option is automatically false in this case. In all other cases, coupled induction is strongly recommended. Only turn it off if the induction energy is not going to be used.
Type: boolean
Default: true
SAPT(DFT)¶
In general, SAPT(DFT) should provide more accurate interaction energy components, and overall interaction energies, than SAPT0. The drawback is SAPT(DFT) method is more computationally demanding than SAPT0, SAPT(DFT) can still be applied to medium-sized or large systems. The SAPT(DFT) module was employed successfully in computations of systems with up to 2000 basis functions, and the code should be scalable to 3000 basis functions. Like higher-order SAPT, SAPT(DFT) requires sufficient memory to hold \(2ovN_aux\) doubles.
SAPT(DFT) requires a few special treatments to obtain accurate
result. The DFT functionals used in SAPT(DFT) need to be asymptotically
corrected with Gradient Regulated Asymptotic Correction scheme (GRAC),
in order to recover the correct long-range asymptotic behavior
(approaching \(-1/r\) as \(r\) approaches infinity). The program
requires manual input of GRAC shift parameter for each monomer through
keywords SAPT_DFT_GRAC_SHIFT_A and SAPT_DFT_GRAC_SHIFT_B,
which should be equal to the difference of the actual ionization
potential and the corresponding Kohn-Sham HOMO energy. The dispersion
term needs to be computed with orbital response for good accuracy,
and it is recommended to enable SAPT_DFT_DO_HYBRID (set to
True by default). The coupled exchange-dispersion energy is usually
estimated by scaling from the uncoupled value either by a fitted fixed
value (suggested initially by [Hesselmann:2014:094107] for a local Hartree–Fock (LHF) formulation and then revised
by [Xie:2022:024801] for non-LHF) or
by the ratio of
coupled and uncoupled dispersion energy (suggested by [Podeszwa:2006:400] ).
This can be controlled by keyword SAPT_DFT_EXCH_DISP_SCALE_SCHEME,
with FIXED using the Hesselmann/Xie approach (PSI4‘s default prior
to Nov 2022), DISP using the Podeszwa approach (PSI4‘s default after Nov 2022),
or NONE for not scaling and using the uncoupled exchange-dispersion
energy directly.
Warning
Since Nov 2022, the defaults of options SAPT_DFT_EXCH_DISP_SCALE_SCHEME and SAPT_DFT_EXCH_DISP_FIXED_SCALE
have been changed. Before, the former defaulted to FIXED with Hesselmann value of 0.686 for the latter. Now, the former defaults to DISP and should you instead select FIXED, the default for the latter is the Xie value of 0.770. This might cause
an older version of PSI4 to produce a different value of
exchange-dispersion energy from the latest version.
Basic Keywords for SAPT(DFT)¶
SAPT_DFT_GRAC_SHIFT_A¶
Monomer A GRAC shift in Hartree
Type: double
Default: 0.0
SAPT_DFT_GRAC_SHIFT_B¶
Monomer B GRAC shift in Hartree
Type: double
Default: 0.0
SAPT_DFT_DO_DHF¶
Compute the Delta-HF correction?
Type: boolean
Default: true
SAPT_DFT_EXCH_DISP_SCALE_SCHEME¶
Scheme for approximating exchange-dispersion for SAPT-DFT. Previous to Nov 2022, default was
FIXEDwith Hesselmann value.NONEUse unscaledExch-Disp2,u.FIXEDUse a fixed factor SAPT_DFT_EXCH_DISP_FIXED_SCALE to scaleExch-Disp2,u.DISPUse the ratio ofDisp2,randDisp2,uto scaleExch-Disp2,u.
Type: string
Possible Values: NONE, FIXED, DISP
Default: FIXED
Advanced Keywords for SAPT(DFT)¶
SAPT_DFT_FUNCTIONAL¶
Underlying funcitonal to use for SAPT(DFT)
Type: string
Default: PBE0
SAPT_DFT_DO_HYBRID¶
Enables the hybrid xc kernel in dispersion?
Type: boolean
Default: true
SAPT_DFT_EXCH_DISP_FIXED_SCALE¶
Exch-disp scaling factor for FIXED scheme for SAPT_DFT_EXCH_DISP_SCALE_SCHEME Default value of 0.770 suggested in Y. Xie, D. G. A. Smith and C. D. Sherrill, 2022 (submitted). Previous to Nov 2022, default value was 0.686 suggested by Hesselmann and Korona, J. Chem. Phys. 141, 094107 (2014).
Type: double
Default: 0.770
SAPT_DFT_MP2_DISP_ALG¶
Which MP2 Exch-Disp module to use?
Type: string
Possible Values: FISAPT, SAPT
Default: SAPT
SAPT_QUIET¶
Interior option to clean up printing
Type: boolean
Default: false
Higher-Order SAPT¶
For smaller systems (up to the size of a nucleic acid base pair), more accurate interaction energies can be obtained through higher-order SAPT computations. The SAPT module can perform density-fitted evaluations of SAPT2, SAPT2+, SAPT2+(3), and SAPT2+3 energies for closed-shell systems only. Publications resulting from the use of the higher-order SAPT code should cite the following: [Hohenstein:2010:014101].
For methods SAPT2+ and above, one can replace the many-body treatment of
dispersion by an improved method based on coupled-cluster doubles (CCD).
This approach tends to give good improvements when dispersion effects
are very large, as in the PCCP dimer (see [Hohenstein:2011:2842]).
As shown in [Parker:2014:094106], whether or not CCD dispersion offers
more accurate interaction energies tends to depend on the SAPT truncation
and basis set employed, due to cancellations of errors. Thanks to
natural orbital methods [Parrish:2013:174102], the SAPT code
is able to include CCD dispersion with only a modest additional cost.
Computations employing CCD dispersion should cite [Parrish:2013:174102].
To request CCD dispersion treatment in a SAPT computation, simply append
(ccd) to the name of the method, as in the following examples
energy('sapt2+(ccd)')
energy('sapt2+(3)(ccd)')
energy('sapt2+3(ccd)')
The \(\delta_{MP2}\) corrections can also be computed automatically
by appending dmp2 to the name of the method, with or without CCD dispersion
energy('sapt2+dmp2')
energy('sapt2+(3)dmp2')
energy('sapt2+3dmp2')
energy('sapt2+(ccd)dmp2')
energy('sapt2+(3)(ccd)dmp2')
energy('sapt2+3(ccd)dmp2')
A brief note on memory usage: the higher-order SAPT code assumes that certain quantities can be held in core. This code requires sufficient memory to hold \(3o^2v^2+v^2N_{aux}\) arrays in core. With this requirement computations on the adenine-thymine complex can be performed with an aug-cc-pVTZ basis in less than 64GB of memory.
Higher-order SAPT is treated separately from the highly optimized SAPT0 code, therefore, higher-order SAPT uses a separate set of keywords. The following keywords are relevant for higher-order SAPT.
Basic Keywords for Higher-order SAPT¶
BASIS¶
Primary basis set, describes the monomer molecular orbitals
Type: string
Possible Values: basis string
Default: No Default
DF_BASIS_SAPT¶
Auxiliary basis set for SAPT density fitting computations. Defaults to a RI basis.
Type: string
Possible Values: basis string
Default: No Default
FREEZE_CORE¶
Specifies how many core orbitals to freeze in correlated computations.
TRUEor1will default to freezing the previous noble gas shell on each atom. In case of positive charges on fragments, an additional shell may be unfrozen, to ensure there are valence electrons in each fragment. WithFALSEor0, no electrons are frozen (with the exception of electrons treated by an ECP). With-1,-2, and-3, the user might request strict freezing of the previous first/second/third noble gas shell on every atom. In this case, when there are no valence electrons, the code raises an exception. More precise control over the number of frozen orbitals can be attained by using the keywords NUM_FROZEN_DOCC (gives the total number of orbitals to freeze, program picks the lowest-energy orbitals) or FROZEN_DOCC (gives the number of orbitals to freeze per irreducible representation) or by the optionPOLICYin combination with appropriate inputs to FREEZE_CORE_POLICY At present,POLICYis an experimental option and is subject to change.
Type: string
Possible Values: FALSE, TRUE, 1, 0, -1, -2, -3, POLICY
Default: FALSE
PRINT¶
The amount of information to print to the output file for the sapt module. For 0, only the header and final results are printed. For 1, (recommended for large calculations) some intermediate quantities are also printed.
Type: integer
Default: 1
Advanced Keywords for Higher-order SAPT¶
DO_CCD_DISP¶
Do CCD dispersion correction in SAPT2+, SAPT2+(3) or SAPT2+3?
Type: boolean
Default: false
DO_MBPT_DISP¶
Do MBPT dispersion correction in SAPT2+, SAPT2+(3) or SAPT2+3, if also doing CCD?
Type: boolean
Default: true
DO_THIRD_ORDER¶
Do compute third-order corrections?
Type: boolean
Default: false
INTS_TOLERANCE¶
Schwarz screening threshold. Minimum absolute value below which all three-index DF integrals and those contributing to four-index integrals are neglected. The default is conservative, but there isn’t much to be gained from loosening it, especially for higher-order SAPT.
Type: conv double
Default: 1.0e-12
SAPT_MEM_CHECK¶
Do force SAPT2 and higher to die if it thinks there isn’t enough memory? Turning this off is ill-advised.
Type: boolean
Default: true
DEBUG¶
The amount of information to print to the output file
Type: integer
Default: 0
MP2 Natural Orbitals¶
One of the unique features of the SAPT module is its ability to use MP2 natural orbitals (NOs) to speed up the evaluation of the triples contribution to dispersion. By transforming to the MP2 NO basis, we can throw away virtual orbitals that are expected to contribute little to the dispersion energy. Speedups in excess of \(50 \times\) are possible. In practice, this approximation is very good and should always be applied. Publications resulting from the use of MP2 NO-based approximations should cite the following: [Hohenstein:2010:104107].
Basic Keywords Controlling MP2 NO Approximations¶
NAT_ORBS_T2¶
Do use MP2 natural orbital approximations for the \(v^4\) block of two-electron integrals in the evaluation of second-order T2 amplitudes? Recommended true for all SAPT computations.
Type: boolean
Default: true
NAT_ORBS_T3¶
Do natural orbitals to speed up evaluation of the triples contribution to dispersion by truncating the virtual orbital space? Recommended true for all SAPT computations.
Type: boolean
Default: true
NAT_ORBS_V4¶
Do use MP2 natural orbital approximations for the \(v^4\) block of two-electron integrals in the evaluation of CCD T2 amplitudes? Recommended true for all SAPT computations.
Type: boolean
Default: true
OCC_TOLERANCE¶
Minimum occupation (eigenvalues of the MP2 OPDM) below which virtual natural orbitals are discarded for in each of the above three truncations
Type: conv double
Default: 1.0e-6
Charge-Transfer in SAPT¶
It is possible to obtain the stabilization energy of a complex due to charge-transfer effects from a SAPT computation. The charge-transfer energy can be computed with the SAPT module as described by Stone and Misquitta [Misquitta:2009:201].
Charge-transfer energies can be obtained from the following calls to the energy function.
energy('sapt0-ct')
energy('sapt2-ct')
energy('sapt2+-ct')
energy('sapt2+(3)-ct')
energy('sapt2+3-ct')
energy('sapt2+(ccd)-ct')
energy('sapt2+(3)(ccd)-ct')
energy('sapt2+3(ccd)-ct')
For now, charge transfer computations are not available with open-shell SAPT0.
A SAPT charge-transfer analysis will perform 5 HF computations: the dimer in the dimer basis, monomer A in the dimer basis, monomer B in the dimer basis, monomer A in the monomer A basis, and monomer B in the monomer B basis. Next, it performs two SAPT computations, one in the dimer basis and one in the monomer basis. Finally, it will print a summary of the charge-transfer results:
SAPT Charge Transfer Analysis
------------------------------------------------------------------------------------------------
SAPT Induction (Dimer Basis) -2.0970 [mEh] -1.3159 [kcal/mol] -5.5057 [kJ/mol]
SAPT Induction (Monomer Basis) -1.1396 [mEh] -0.7151 [kcal/mol] -2.9920 [kJ/mol]
SAPT Charge Transfer -0.9574 [mEh] -0.6008 [kcal/mol] -2.5137 [kJ/mol]
These results are for the water dimer geometry shown above computed with SAPT0/aug-cc-pVDZ.
Monomer-Centered Basis Computations¶
The charge-transfer analysis above is carried out by taking the
difference between SAPT induction as calculated in the dimer-centered
basis (i.e., each monomer sees the basis functions on both monomers)
vs. the monomer-centered basis (i.e., each monomer utilizes only its
own basis set). It is also possible to run a closed-shell SAPT computation at any
level using only the monomer-centered basis. To do this, simply add
sapt_basis='monomer' to the energy function, such as
energy('sapt2',sapt_basis='monomer')
This procedure leads to faster compuations, but it converges more slowly towards the complete basis set limit than the default procedure, which uses the dimer-centered basis set. Hence, monomer-centered basis SAPT computations are not recommended. The open-shell SAPT0 code is not compatible yet with monomer-centered computations.
Computations with Mid-bonds¶
SAPT computations with midbonds can be accomplished by adding a third ghost monomer to the computation. For example
molecule dimer {
0 1
He 0 0 5
--
0 1
He 0 0 -5
--
0 1
@He 0 0 0
}
Here the functions of the third monomer will be added to the virtual space of the entire computation. Note that an error will be thrown if each atom in the third monomer is not a ghost to prevent confusion with three-body SAPT which is not currently supported by Psi4.
Interpreting SAPT Results¶
We will examine the results of a SAPT2+3/aug-cc-pVDZ computation on the water dimer. This computation can be performed with the following input:
molecule water_dimer {
0 1
O -1.551007 -0.114520 0.000000
H -1.934259 0.762503 0.000000
H -0.599677 0.040712 0.000000
--
0 1
O 1.350625 0.111469 0.000000
H 1.680398 -0.373741 -0.758561
H 1.680398 -0.373741 0.758561
units angstrom
}
set globals {
basis aug-cc-pvdz
guess sad
scf_type df
}
set sapt {
print 1
nat_orbs_t2 true
freeze_core true
}
energy('sapt2+3')
To reiterate some of the options mentioned above: the NAT_ORBS_T2 option will compute MP2 natural orbitals and use them in the evaluation of the triples correction to dispersion, and the FREEZE_CORE option will freeze the core throughout the SAPT computation. This SAPT2+3/aug-cc-pVDZ computation produces the following results:
SAPT Results
--------------------------------------------------------------------------------------------------------
Electrostatics -13.06509118 [mEh] -8.19846883 [kcal/mol] -34.30239689 [kJ/mol]
Elst10,r -13.37542977 [mEh] -8.39320925 [kcal/mol] -35.11719087 [kJ/mol]
Elst12,r 0.04490350 [mEh] 0.02817737 [kcal/mol] 0.11789413 [kJ/mol]
Elst13,r 0.26543510 [mEh] 0.16656305 [kcal/mol] 0.69689985 [kJ/mol]
Exchange 13.41768202 [mEh] 8.41972294 [kcal/mol] 35.22812415 [kJ/mol]
Exch10 11.21822294 [mEh] 7.03954147 [kcal/mol] 29.45344432 [kJ/mol]
Exch10(S^2) 11.13802706 [mEh] 6.98921779 [kcal/mol] 29.24289005 [kJ/mol]
Exch11(S^2) 0.04558907 [mEh] 0.02860757 [kcal/mol] 0.11969410 [kJ/mol]
Exch12(S^2) 2.15387002 [mEh] 1.35157390 [kcal/mol] 5.65498573 [kJ/mol]
Induction -3.91313050 [mEh] -2.45552656 [kcal/mol] -10.27392413 [kJ/mol]
Ind20,r -4.57530818 [mEh] -2.87104935 [kcal/mol] -12.01247162 [kJ/mol]
Ind30,r -4.91714746 [mEh] -3.08555675 [kcal/mol] -12.90997067 [kJ/mol]
Ind22 -0.83718642 [mEh] -0.52534243 [kcal/mol] -2.19803293 [kJ/mol]
Exch-Ind20,r 2.47828501 [mEh] 1.55514739 [kcal/mol] 6.50673730 [kJ/mol]
Exch-Ind30,r 4.33916119 [mEh] 2.72286487 [kcal/mol] 11.39246770 [kJ/mol]
Exch-Ind22 0.45347471 [mEh] 0.28455969 [kcal/mol] 1.19059785 [kJ/mol]
delta HF,r (2) -1.43239563 [mEh] -0.89884187 [kcal/mol] -3.76075473 [kJ/mol]
delta HF,r (3) -0.85440936 [mEh] -0.53614999 [kcal/mol] -2.24325177 [kJ/mol]
Dispersion -3.62000698 [mEh] -2.27158877 [kcal/mol] -9.50432831 [kJ/mol]
Disp20 -3.54291925 [mEh] -2.22321549 [kcal/mol] -9.30193450 [kJ/mol]
Disp30 0.05959979 [mEh] 0.03739944 [kcal/mol] 0.15647926 [kJ/mol]
Disp21 0.11216169 [mEh] 0.07038252 [kcal/mol] 0.29448051 [kJ/mol]
Disp22 (SDQ) -0.17892163 [mEh] -0.11227502 [kcal/mol] -0.46975875 [kJ/mol]
Disp22 (T) -0.47692534 [mEh] -0.29927518 [kcal/mol] -1.25216749 [kJ/mol]
Est. Disp22 (T) -0.54385233 [mEh] -0.34127251 [kcal/mol] -1.42788430 [kJ/mol]
Exch-Disp20 0.64545587 [mEh] 0.40502969 [kcal/mol] 1.69464439 [kJ/mol]
Exch-Disp30 -0.01823410 [mEh] -0.01144207 [kcal/mol] -0.04787362 [kJ/mol]
Ind-Disp30 -0.91816882 [mEh] -0.57615966 [kcal/mol] -2.41065224 [kJ/mol]
Exch-Ind-Disp30 0.76487181 [mEh] 0.47996433 [kcal/mol] 2.00817094 [kJ/mol]
Total HF -5.68662563 [mEh] -3.56841161 [kcal/mol] -14.93023559 [kJ/mol]
Total SAPT0 -8.58408901 [mEh] -5.38659740 [kcal/mol] -22.53752571 [kJ/mol]
Total SAPT2 -6.72343814 [mEh] -4.21902130 [kcal/mol] -17.65238683 [kJ/mol]
Total SAPT2+ -7.33405042 [mEh] -4.60218631 [kcal/mol] -19.25554938 [kJ/mol]
Total SAPT2+(3) -7.00901553 [mEh] -4.39822383 [kcal/mol] -18.40217026 [kJ/mol]
Total SAPT2+3 -7.18054663 [mEh] -4.50586123 [kcal/mol] -18.85252518 [kJ/mol]
Special recipe for scaled SAPT0 (see Manual):
Electrostatics sSAPT0 -13.37542977 [mEh] -8.39320925 [kcal/mol] -35.11719087 [kJ/mol]
Exchange sSAPT0 11.21822294 [mEh] 7.03954147 [kcal/mol] 29.45344432 [kJ/mol]
Induction sSAPT0 -3.47550008 [mEh] -2.18090932 [kcal/mol] -9.12492546 [kJ/mol]
Dispersion sSAPT0 -2.88342055 [mEh] -1.80937379 [kcal/mol] -7.57042064 [kJ/mol]
Total sSAPT0 -8.51612746 [mEh] -5.34395089 [kcal/mol] -22.35909265 [kJ/mol]
--------------------------------------------------------------------------------------------------------
At the bottom of this output are the total SAPT energies (defined above),
they are composed of subsets of the individual terms printed above. The
individual terms are grouped according to the component of the interaction
to which they contribute. The total component energies (i.e.,
electrostatics, exchange, induction, and dispersion) represent what we
regard as the best estimate available at a given level of SAPT computed
from a subset of the terms of that grouping. The groupings shown above are
not unique and are certainly not rigorously defined. We regard the groupings
used in PSI4 as a “chemist’s grouping” as opposed to a more
mathematically based grouping, which would group all exchange terms
(i.e. \(E_{exch-ind,resp}^{(20)}\), \(E_{exch-disp}^{(20)}\), etc.) in
the exchange component. A final note is that both Disp22(T)
and Est.Disp22(T) results appear if MP2 natural orbitals are
used to evaluate the triples correction to dispersion. The Disp22(T)
result is the triples correction as computed in the truncated NO basis;
Est.Disp22(T) is a scaled result that attempts to recover
the effect of the truncated virtual space and is our best estimate. The Est.Disp22(T)
value is used in the SAPT energy and dispersion component (see [Hohenstein:2010:104107]
for details). Finally, this part of the output file contains sSAPT0, a special scaling
scheme of the SAPT0 energy that can yield improved results and was described in more details
above. The corresponding scaled total component energies are printed as well.
As mentioned above, SAPT results with scaled exchange are also optionally available by setting the EXCH_SCALE_ALPHA keyword. When activated, the unscaled results are printed first as reported above, and then repeated with exchange scaling for all relevant terms:
SAPT Results ==> ALL S2 TERMS SCALED (see Manual) <==
Scaling factor (Exch10/Exch10(S^2))^{Alpha} = 1.007200
with Alpha = 1.000000
--------------------------------------------------------------------------------------------------------
Electrostatics -13.06509118 [mEh] -8.19846883 [kcal/mol] -34.30239689 [kJ/mol]
Elst10,r -13.37542977 [mEh] -8.39320925 [kcal/mol] -35.11719087 [kJ/mol]
Elst12,r 0.04490350 [mEh] 0.02817737 [kcal/mol] 0.11789413 [kJ/mol]
Elst13,r 0.26543510 [mEh] 0.16656305 [kcal/mol] 0.69689985 [kJ/mol]
Exchange sc. 13.43351854 [mEh] 8.42966050 [kcal/mol] 35.26970292 [kJ/mol]
Exch10 11.21822294 [mEh] 7.03954147 [kcal/mol] 29.45344432 [kJ/mol]
Exch10(S^2) 11.13802706 [mEh] 6.98921779 [kcal/mol] 29.24289005 [kJ/mol]
Exch11(S^2) sc. 0.04591732 [mEh] 0.02881355 [kcal/mol] 0.12055592 [kJ/mol]
Exch12(S^2) sc. 2.16937828 [mEh] 1.36130548 [kcal/mol] 5.69570268 [kJ/mol]
Induction sc. -3.90986540 [mEh] -2.45347768 [kcal/mol] -10.26535160 [kJ/mol]
Ind20,r -4.57530818 [mEh] -2.87104935 [kcal/mol] -12.01247162 [kJ/mol]
Ind30,r -4.91714746 [mEh] -3.08555675 [kcal/mol] -12.90997067 [kJ/mol]
Ind22 -0.83718642 [mEh] -0.52534243 [kcal/mol] -2.19803293 [kJ/mol]
Exch-Ind20,r sc. 2.49612913 [mEh] 1.56634474 [kcal/mol] 6.55358703 [kJ/mol]
Exch-Ind30,r sc. 4.37040396 [mEh] 2.74247000 [kcal/mol] 11.47449560 [kJ/mol]
Exch-Ind22 sc. 0.45673981 [mEh] 0.28660857 [kcal/mol] 1.19917038 [kJ/mol]
delta HF,r (2) sc. -1.45023975 [mEh] -0.91003922 [kcal/mol] -3.80760445 [kJ/mol]
delta HF,r (3) sc. -0.90349624 [mEh] -0.56695248 [kcal/mol] -2.37212939 [kJ/mol]
Dispersion sc. -3.60998364 [mEh] -2.26529903 [kcal/mol] -9.47801205 [kJ/mol]
Disp20 -3.54291925 [mEh] -2.22321549 [kcal/mol] -9.30193450 [kJ/mol]
Disp30 0.05959979 [mEh] 0.03739944 [kcal/mol] 0.15647926 [kJ/mol]
Disp21 0.11216169 [mEh] 0.07038252 [kcal/mol] 0.29448051 [kJ/mol]
Disp22 (SDQ) -0.17892163 [mEh] -0.11227502 [kcal/mol] -0.46975875 [kJ/mol]
Disp22 (T) -0.47692534 [mEh] -0.29927518 [kcal/mol] -1.25216749 [kJ/mol]
Est. Disp22 (T) -0.54385233 [mEh] -0.34127251 [kcal/mol] -1.42788430 [kJ/mol]
Exch-Disp20 sc. 0.65010327 [mEh] 0.40794598 [kcal/mol] 1.70684615 [kJ/mol]
Exch-Disp30 sc. -0.01836538 [mEh] -0.01152445 [kcal/mol] -0.04821832 [kJ/mol]
Ind-Disp30 -0.91816882 [mEh] -0.57615966 [kcal/mol] -2.41065224 [kJ/mol]
Exch-Ind-Disp30 sc. 0.77037903 [mEh] 0.48342016 [kcal/mol] 2.02263015 [kJ/mol]
Total HF -5.68662563 [mEh] -3.56841161 [kcal/mol] -14.93023559 [kJ/mol]
Total SAPT0 sc. -8.57944161 [mEh] -5.38368112 [kcal/mol] -22.52532395 [kJ/mol]
Total SAPT2 sc. -6.69968912 [mEh] -4.20411857 [kcal/mol] -17.59003378 [kJ/mol]
Total SAPT2+ sc. -7.31030140 [mEh] -4.58728357 [kcal/mol] -19.19319632 [kJ/mol]
Total SAPT2+(3) sc. -6.98526650 [mEh] -4.38332109 [kcal/mol] -18.33981720 [kJ/mol]
Total SAPT2+3 sc. -7.15142168 [mEh] -4.48758504 [kcal/mol] -18.77605762 [kJ/mol]
--------------------------------------------------------------------------------------------------------
The scaling factor is reported at the top (here 1.0072) together with the
\(\alpha\) parameter. All terms that are scaled are indicated by the sc.
label. Note that if Exch10 is less than \(10^{-5}\), the scaling factor is
set to \(1.0\).
Caution
To density fit the dispersion terms in SAPT, the RI auxiliary basis set (e.g., aug-cc-pVDZ-RI) controlled through DF_BASIS_SAPT performs well. For Fock-type terms (i.e., electrostatics, exchange, induction, and core Fock matrix elements in exchange-dispersion), the JKFIT density-fitting auxiliary basis (e.g., aug-cc-pVDZ-JKFIT) is more appropriate. The FISAPT module has always used JKFIT in this role. The SAPT module newly (see fitting notes ) uses JKFIT for computations targeting SAPT0 and sSAPT0 methods. But the SAPT module still uses the RI basis for higher-order SAPT. For heavier elements (i.e., second-row and beyond), the RI auxiliary basis is unsound for this role (insufficiently flexible). For higher-order methods in SAPT module, there is no workaround; on-the-fly construction of an auxiliary basis through Cholesky decomposition (not implemented) is the long-term solution.
Spin-Flip SAPT¶
SAPT0 with two open-shell references will always yield a high-spin complex. In order to obtain a SAPT-based estimate of the splittings between different spin states of a complex the first-order exchange energies for all multiplets can be shown to be a linear combination of two matrix elements: a diagonal exchange term that determines the spin-averaged effect and a spin-flip term responsible for the splittings between the states. The numerical factors in this linear combination are determined solely by the Clebsch-Gordan coefficients: accordingly, the \(S^{2}\) approximation implies a Heisenberg Hamiltonian picture with a single coupling strength parameter determining all the splittings. This method can be invoked with energy(“SF-SAPT”) and publications resulting from the use of the SF-SAPT code should cite the following publications: [Patkowski:2018:164110]
Higher-Order Exchange Terms without Single-Exchange Approximation¶
Recently, several SAPT higher-order exchange terms have been derived without the \(S^{2}\) approximation: \(E_{exch-ind}^{(20)}\) [Schaffer:2012:1235], \(E_{exch-disp}^{(20)}\) [Schaffer:2013:2570], and \(E_{exch-ind}^{(30)}\) [Waldrop:2021:024103]. The second-order terms can be computed with the following settings:
set SAPT_DFT_FUNCTIONAL HF
set DO_IND_EXCH_SINF true # calculate Exch-Ind20 (S^inf)
set SAPT_DFT_MP2_DISP_ALG fisapt
set DO_DISP_EXCH_SINF true # calculate Exch-Disp20 (S^inf)
energy('sapt(dft)')
and the third-order exchange-induction term is computed as follows:
set DO_IND30_EXCH_SINF true # calculate Exch-Ind30 (S^inf)
energy('sapt2+3')
These calculations are performed with the atomic orbital and density-fitting scheme described in the Supplementary Material to [Smith:2020:184108] for the second-order terms and in [Waldrop:2021:024103] for the third-order exchange induction. The coupled (response) version of the exchange-induction corrections are also calculated, exactly for \(E_{exch-ind,resp}^{(20)}\) and by scaling the uncoupled term for \(E_{exch-ind,resp}^{(30)}\).
S^inf Keywords¶
DO_IND_EXCH_SINF¶
For SAPT0 or SAPT(DFT), compute the non-approximated second-order exchange-induction term.
Type: boolean
Default: false
DO_DISP_EXCH_SINF¶
For SAPT0 or SAPT(DFT), compute the non-approximated second-order exchange-dispersion term.
Type: boolean
Default: false
DO_IND30_EXCH_SINF¶
For SAPT2+3, compute the non-approximated third-order exchange-induction term.
Type: boolean
Default: false
SAPT0-D¶
In SAPT0, the computation of \(E_{disp}^{(20)} + E_{exch-disp}^{(20)}\) represents the computational bottleneck. One can avoid this bottleneck by replacing these dispersion terms with the empirical D3 corrections developed by Grimme.
Grimme’s dispersion corrections are discussed here.
The corresponding method, termed SAPT0-D, thus relies on empirically fit parameters specific to SAPT0/jun-cc-pVDZ. While SAPT0-D can be used with any of the -D variants using default parameters optimized for Hartee–Fock interaction energies, we recommend using the refit parameters with Becke-Johnson damping, as described in [Schriber:2021:234107]. Again, use of SAPT0-D with a basis set other than jun-cc-pVDZ is not tested and not guaranteed to give meaningful results without refitting the dispersion parameters. A simple water dimer computation using SAPT0-D may look like:
molecule water_dimer {
0 1
O -1.551007 -0.114520 0.000000
H -1.934259 0.762503 0.000000
H -0.599677 0.040712 0.000000
--
0 1
O 1.350625 0.111469 0.000000
H 1.680398 -0.373741 -0.758561
H 1.680398 -0.373741 0.758561
units angstrom
no_reorient
symmetry c1
}
set basis jun-cc-pvdz
energy('sapt0-d3mbj') # runs the recommended dispersion correction
energy('sapt0-d3') # tests an alternative damping scheme/parameterization
Given the naturally pairwise-atomic nature of these empirical dispersion corrections,
integration with existing FSAPT functionality is also available simply by calling
energy("fsapt0-d3mbj"). See FSAPT documentation for more details on using FSAPT
for functional group analyses.
