Searching....

# 12.14.4 Running an XSAPT+MBD Job

(July 10, 2023)

In order to run an XSAPT job, at least two separate jobs must be run, a third is optional and suggested in the case where the induction energy in the system is supsected of being important. The first job that must be run is a second-order SAPT job without charge embedding on the full system. This job provides the value for the electrostatic term, the exchange term, and the dispersion term. (For an XSAPT+MBD method 174 Carter-Fenk K. et al.
J. Phys. Chem. Lett.
(2019), 10, pp. 2706.
job, the dispersion value will be the dispersion value obtained from the MBD calculation). The second job that needs to be run is a SAPT job with charge embedding on the full system. This SAPT job with charge embedding acts as a correction for the induction energy and must be run on each pair of monomers. The final induction energy for an XSAPT job is then the induction energy from the SAPT job without charge embedding plus the difference between the total SAPT energy with charge embedding minus the total SAPT energy without charge embedding. In equation form, the XSAPT contribution to the induction energy would be:

 $E_{pol}=\sum_{A}\sum_{B>A}\left[E_{AB}^{\text{XSAPT}}(AB)-E_{AB}^{\text{SAPT}}\right]$ (12.66)

Here the term in brackets is the difference between the charge-embedded energy for the $AB$ dimer and the SAPT energy computed without charge embedding. Therefore, the electrostatic, exchange, and dispersion energy would be the equivalent to the same terms calculated from a normal SAPT calculation, it is the induction energy that changes upon an XSAPT calculation. An additional calculation can be performed to further correct the induction energy, this is the $\delta$HF calculation, which accounts for higher than second-order induction energy. The resulting delta(SCF) energy term is then added to the induction energy term calculated as described above.

The above procedure is the required procedure for performing an XSAPT ($+\delta$HF) calculaiton on a dimer complex. In order to perform an XAPT ($+\delta$HF) calculation on an n-mer complex, several more calculations need to be performed. The first of which is to run a SAPT job with charge embedding on each pair of monomers. This is required because the polarization of the dimers differs from that of the dimers embedded in the charges of the entire system. The second set of jobs, if the $\delta$HF correction is being performed, is a $\delta$HF job on every pair of monomers. These sets of jobs will again provide a correction to the induction energy for the system.

Two scripts to help the user perform an XSAPT+MBD+$\delta$HF calculation exist in the $QC/bin/tools directory. The script named xsapt_creation.py creates three files for an arbitrary geometry (defined by the user), a SAPT job without embedding, a SAPT job with embedding, and a $\delta$HF job. See the script for instructions on how to use and modify the script for personal use. The other script named xsapt_data_collection.py looks through the output files of the three jobs created by xsapt_creation.py and performs the needed arithmetic to obtain the correct energy decomposition for an XSAPT+MBD+$\delta$HF job. See the script for instructions on how to use and modify it for personal use. Example 12.42 AO-XSAPT(KS)+D3 calculation of water-water interaction. $molecule
0 1
--
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
$end$rem
JOBTYPE            xsapt
EXCHANGE           gen
BASIS              aug-cc-pVTZ
MEM_TOTAL          46000
MEM_STATIC         4000
AO2MO_DISK         35000
CHELPG_DX          5
CHELPG_H           110
CHELPG_HA          590
$end$xpol
embed   charges
charges CHELPG
DFT-LRC
print   3
$end$sapt
algorithm     ao  ! for use with +aiD dispersion
order         2   ! 2nd-order SAPT
basis         projected
Dispersion    aiD3
print         3
$end$xc_functional
x   wPBE  1.0
c   PBE   1.0
$end$lrc_omega
502
502
\$end