SAPT0 results unphysically different for isomers?

I recently performed SAPT0 calculations on two nucleic acid base pairs. The base pairs are constitutional isomers. Given the same number of atoms in both base pairs, and given the same number of hydrogen bonds, I’m concerned about the SAPT0 result values. Specifically although the Total SAPT0 seems appropriate when comparing the two base pairs, the absolute values of the individual SAPT terms appear unphysical between the two – e.g., exchange values of 16.4 vs 45.7 kcal/mol, etc. Sooner or later I have to pin the SAPT0 output on molecular features. Perhaps there is something wonky with the input file parameters shown below? The HF results are fine when the two base pairs are compared.


      Psi4: An Open-Source Ab Initio Electronic Structure Package
                           Psi4 1.1 release

                     Git: Rev {HEAD} add49b9 

R. M. Parrish, L. A. Burns, D. G. A. Smith, A. C. Simmonett,
A. E. DePrince III, E. G. Hohenstein, U. Bozkaya, A. Yu. Sokolov,
R. Di Remigio, R. M. Richard, J. F. Gonthier, A. M. James,
H. R. McAlexander, A. Kumar, M. Saitow, X. Wang, B. P. Pritchard,
P. Verma, H. F. Schaefer III, K. Patkowski, R. A. King, E. F. Valeev,
F. A. Evangelista, J. M. Turney, T. D. Crawford, and C. D. Sherrill,
J. Chem. Theory Comput. in press (2017).
(doi: 10.1021/acs.jctc.7b00174)



molecule dimer {
0 1
–snip nucleobase A coords–

0 1
–snip nucleobase B coords–

 units angstrom
 symmetry c1


set globals {
basis aug-cc-pvdz
df_basis_scf aug-cc-pvdz-jkfit
df_basis_sapt aug-cc-pvdz-ri
guess sad
scf_type df



dna base pair 1 results:

SAPT Results 

Electrostatics                -53.83325110 [mEh]     -33.78087648 [kcal/mol]    -141.33920076 [kJ/mol]
  Elst10,r                    -53.83325110 [mEh]     -33.78087648 [kcal/mol]    -141.33920076 [kJ/mol]

Exchange                       26.20929212 [mEh]      16.44657979 [kcal/mol]      68.81249645 [kJ/mol]
  Exch10                       26.20929212 [mEh]      16.44657979 [kcal/mol]      68.81249645 [kJ/mol]
  Exch10(S^2)                  26.04854472 [mEh]      16.34570927 [kcal/mol]      68.39045417 [kJ/mol]

Induction                     -18.28986948 [mEh]     -11.47706685 [kcal/mol]     -48.02005232 [kJ/mol]
  Ind20,r                     -19.21819951 [mEh]     -12.05960276 [kcal/mol]     -50.45738280 [kJ/mol]
  Exch-Ind20,r                  6.87859136 [mEh]       4.31638142 [kcal/mol]      18.05974161 [kJ/mol]
  delta HF,r (2)               -5.95026133 [mEh]      -3.73384551 [kcal/mol]     -15.62241112 [kJ/mol]

Dispersion                    -12.02951381 [mEh]      -7.54863420 [kcal/mol]     -31.58348851 [kJ/mol]
  Disp20                      -14.05565170 [mEh]      -8.82005497 [kcal/mol]     -36.90311354 [kJ/mol]
  Exch-Disp20                   2.02613789 [mEh]       1.27142077 [kcal/mol]       5.31962503 [kJ/mol]
  Disp20 (SS)                  -7.02782585 [mEh]      -4.41002749 [kcal/mol]     -18.45155677 [kJ/mol]
  Disp20 (OS)                  -7.02782585 [mEh]      -4.41002749 [kcal/mol]     -18.45155677 [kJ/mol]
  Exch-Disp20 (SS)              1.25350142 [mEh]       0.78658405 [kcal/mol]       3.29106798 [kJ/mol]
  Exch-Disp20 (OS)              0.77263647 [mEh]       0.48483672 [kcal/mol]       2.02855705 [kJ/mol]

Total HF -45.91382846 [mEh] -28.81136354 [kcal/mol] -120.54675663 [kJ/mol]
Total SAPT0 -57.94334228 [mEh] -36.35999774 [kcal/mol] -152.13024514 [kJ/mol]

Special recipe for scaled SAPT0 (see Manual):
Electrostatics sSAPT0 -53.83325110 [mEh] -33.78087648 [kcal/mol] -141.33920076 [kJ/mol]
Exchange sSAPT0 26.20929212 [mEh] 16.44657979 [kcal/mol] 68.81249645 [kJ/mol]
Induction sSAPT0 -18.16173720 [mEh] -11.39666263 [kcal/mol] -47.68364102 [kJ/mol]
Dispersion sSAPT0 -11.99177154 [mEh] -7.52495056 [kcal/mol] -31.48439618 [kJ/mol]
Total sSAPT0 -57.77746773 [mEh] -36.25590988 [kcal/mol] -151.69474151 [kJ/mol]


dna base pair 2 results:

SAPT Results

Electrostatics                -78.13711547 [mEh]     -49.03178226 [kcal/mol]    -205.14899667 [kJ/mol]
  Elst10,r                    -78.13711547 [mEh]     -49.03178226 [kcal/mol]    -205.14899667 [kJ/mol]

Exchange                       72.83149470 [mEh]      45.70245482 [kcal/mol]     191.21908934 [kJ/mol]
  Exch10                       72.83149470 [mEh]      45.70245482 [kcal/mol]     191.21908934 [kJ/mol]
  Exch10(S^2)                  71.64474227 [mEh]      44.95775640 [kcal/mol]     188.10327082 [kJ/mol]

Induction                     -35.89817406 [mEh]     -22.52644526 [kcal/mol]     -94.25065600 [kJ/mol]
  Ind20,r                     -38.55754529 [mEh]     -24.19522597 [kcal/mol]    -101.23283517 [kJ/mol]
  Exch-Ind20,r                 18.12567024 [mEh]      11.37403027 [kcal/mol]      47.58894721 [kJ/mol]
  delta HF,r (2)              -15.46629900 [mEh]      -9.70524956 [kcal/mol]     -40.60676804 [kJ/mol]

Dispersion                    -20.44753075 [mEh]     -12.83101979 [kcal/mol]     -53.68499197 [kJ/mol]
  Disp20                      -25.08877127 [mEh]     -15.74344231 [kcal/mol]     -65.87056896 [kJ/mol]
  Exch-Disp20                   4.64124052 [mEh]       2.91242252 [kcal/mol]      12.18557699 [kJ/mol]
  Disp20 (SS)                 -12.54438563 [mEh]      -7.87172116 [kcal/mol]     -32.93528448 [kJ/mol]
  Disp20 (OS)                 -12.54438563 [mEh]      -7.87172116 [kcal/mol]     -32.93528448 [kJ/mol]
  Exch-Disp20 (SS)              2.83563971 [mEh]       1.77939086 [kcal/mol]       7.44497207 [kJ/mol]
  Exch-Disp20 (OS)              1.80560081 [mEh]       1.13303166 [kcal/mol]       4.74060492 [kJ/mol]

Total HF -41.20379483 [mEh] -25.85577269 [kcal/mol] -108.18056333 [kJ/mol]
Total SAPT0 -61.65132558 [mEh] -38.68679249 [kcal/mol] -161.86555530 [kJ/mol]

Special recipe for scaled SAPT0 (see Manual):
Electrostatics sSAPT0 -78.13711547 [mEh] -49.03178226 [kcal/mol] -205.14899667 [kJ/mol]
Exchange sSAPT0 72.83149470 [mEh] 45.70245482 [kcal/mol] 191.21908934 [kJ/mol]
Induction sSAPT0 -34.98244899 [mEh] -21.95181907 [kcal/mol] -91.84641982 [kJ/mol]
Dispersion sSAPT0 -20.21305112 [mEh] -12.68388160 [kcal/mol] -53.06936573 [kJ/mol]
Total sSAPT0 -60.50112088 [mEh] -37.96502812 [kcal/mol] -158.84569288 [kJ/mol]

The total IE is only 4 kcal/mol different btwn the isomers, which is reasonable. If isomer 2 has, by its arrangement of fragments A & B, greater steric hinderance and more opportunity for electrostatic attraction, then your results don’t seem shocking.

Input file looks good, provided you don’t have post-Ne elements in there (admittedly unlikely for a nucleobase).

As I mentioned the total energies are fine, and therefore anything derived from them like the IE.

I’m questioning the individual SAPT energies leading up to that (electrostatics, exchange, induction, dispersion). I can get the IE’s from other programs. My focus with psi4 was to learn about how the IE is decomposed into electrostatic, exchange, induction, and dispersion energy components that I can subsequently interpret in molecular terms.

The specific unphysical example I gave above had to do with comparing the Exchange term between the two base pairs. How can the exchange term with one base pair be three times less than the other base pair: 16 kcal/mol vs 45 kcal/mol? All the other terms have equally radically different absolute values as well. This is what I mean by unphysical. Keep in mind that base pairing is edge-edge interactions of flat purines and pyrimidines.

Does one simply not worry about the absolute values of the SAPT electrostatic, exchange, induction, and dispersion terms when making comparisons? If so, I don’t think I can use them in a paper.

Yes, I meant to convey that the large differences in SAPT components between isomers are not immediately alarming, esp. if they reflect your qualitative expectations. At the SAPT0 level, I wouldn’t defend to the death the absolute value of any term (actually, benchmarking SAPT components even at SAPT2+3 suffers from lack of high-level reference calculations), but you can certainly pay attention to trends of the terms.

Perhaps this suppmat showing how SAPT components vary (even 10s of kcal/mol) in base pairs being shifted, twisted, rolled, etc. will be orienting? Relevant tables start at S26.

Looking at the potential plots throughout the supplementary information was helpful, especially those associated with rise. This reinforces my initial intuition to try and interpret the differences in exchange values for base pairing in terms of hydrogen bonding distances.