I try to use the I-SAPT method in order to study from which parts of two molecular units the respective SAPT energy components arise. When I do the plotting I obtain the “difference-wavefunction” before and after perfomring SAPT. In the tutorial with the prototypical diol structure it is said, that electron density is relocated. So which energy component of SAPT does this plot visualize as relocated electron density? Only induction? Or are all components included? Also dispersion? Because dispersion is actually the crucial component I want to prove and study.
You raise a somewhat thorny question, but I’ll try to answer as clearly as possible.
What is being visualized in those cube files is the change in electron density between the zeroth-order monomer wavefunction(s) and the fully interacting monomer wavefunction(s). This will include higher-order induction effects, but not dispersion.
F-/ISAPT performs several SCF cycles, each of which contributes differently to the overall computation:
Isolated monomer A wavefunction
Isolated monomer B wavefunction
Isolated monomer C wavefunction
Monomer A orbital relaxation in presence of Monomer C, constrained orthogonal to Monomer B
Forms zeroth-order monomer A wavefunction
Monomer B orbital relaxation in presence of Monomer C, constrained orthogonal to Monomer A
Forms zeroth-order monomer B wavefunction
Fully interacting supersystem wavefunction (monomers A, B, & C all together)
Used to compute the dHF correction
The difference in electron density between the zeroth-order wavefunctions for monomer A and B (formed in SCF cycles #3 & #4 above) and their density once localized from the fully interacting wavefunction (formed in SCF #6 above) is what is plotted in the cube files. Because all of this embedding is done at the Hartree–Fock level, it will include electrostatics, exchange, and higher-order induction, but it will unfortunately not include any dispersion.
I am pretty certain that there isn’t currently any mechanism available in the literature to visualize differences in electron density due exclusively to dispersion. You might be able to look at differences in natural orbital populations between the zeroth-order wavefunctions and fully interacting ones for a “dispersion-aware” change in electron density, but I’m not even sure if that is completely kosher.
If you believe DFT results, you can run any of the VV10 dispersion functionals and do the VV10-part post-SCF and self-consistently. The difference in their density is then due to dispersion. To my knowledge these density differences are tiny for equilibrium distances.
thank you so much for the extensive explanation! Now I know a lot more!
However, another question came into my mind:
When I perform F-SAPT with two non-covalently connected, but interacting molecules, then I do get a visualization file for the dispersion component or not? It is shown in this F-SAPT tutorial, that was made Rob Parrish I think. Why is this then not possible for I-SAPT? Or do I get something wrong?
Inside the fsapt/ subdirectory output by Psi4 (for either F-SAPT or ISAPT computations), you do get a Disp.dat file. Feeding that through the PyMOL visualization scripts will show you the dispersion interaction energies between the indicated fragments, but that’s different than visualizing the difference in electron density due to the dispersion interactions.
I suppose the answer comes down to the question: which are you hoping to visualize? The interaction energies from the different SAPT components or the difference wavefunction itself?
Hm, that’s a really interesting point. @jeffschriber and I were inspired by this question and have been thinking about how we might be able to build a “density” from the dispersion amplitudes themselves, but we completely forgot to consider VV10. Thanks, @hokru!
I come from a kind of physical organic chemistry direction and computations are rather a tool for me in order to predict or explain certain molecular properties. For example, if a certain conformer of a dimer species is stabilized by dispersion and from which moieties this interaction come from.
So I think what can be done with the Disp.dat is what I am looking for, but can I use that also for I-SAPT, when my dimer is an intramolecular one and has a linker C?
No worries, there’s a great number of technical details in some of these methods and software that can be confusing, even for folks that develop this software.
You can absolutely get each of the Elst.dat, Exch.dat, IndAB.dat, IndBA.dat, and Disp.dat files in the fsapt/ subdirectory from an ISAPT computation. This is because both F-SAPT and ISAPT are computed by the same code, so it will produce the same output regardless of whether there’s a linker.
The only thing you’ll have to keep in mind is that when you are visualizing the dispersion energies (or any other IE component using the PyMOL visualization scripts and .dat files) you won’t be able to see the interactions between any moieties from monomer A or B and the linker, C. You’ll only see the interactions
between moieties on A and B.
Hi everyone, thanks a lot for you comments. I know that this post is old but I have more questions if anybody can help me please.
Do you know what’s are the interactions energy that I can observe in the directory /psi4conda/plot? The prefix of the .cube are named V_.cube, dH_.cube and D_.cube. The underscore is/are the fragments. Is there an example in VMD to get more insight? I did a visualization and I am pretty sure what is the electrostatic component visualization but the other cubes I don’t have any idea.