Question about SAPT's energy contributions

Dear colleagues, I am a new member of the group.

I have a question about SAPT.

I varied the distance between the carbon atoms of two methane molecules and obtained the energies via SAPT.

I know that with the exchange + dispersion energy contributions I can calculate the Lennard-Jones interaction parameters by adjusting the curve obtained in psi4 to the model, for example ( lennard-jones-potential-with-quantum-mechanical-methods).

Regarding the eletrostatics + induction energy curve, could I use the values obtained to calculate the charge of the carbon atom in some way?

My idea was to try to fit the curve to Coulomb’s potential. I honestly don’t know if that makes sense.

I don’t know if this kind of doubt is allowed on the forum. If not, I’m sorry.

Thanks in advance.

I think there are two problems with the idea of trying to fit a partial charge to the SAPT electrostatic+induction curve. First, the SAPT developers have indicated that the accuracy of the SAPT components is likely lower than the total interaction energy (see this post and this post). This is because the accuracy of SAPT relies on error calcellation between the various terms in the SAPT decomposition, and this is especially true for SAPT0.

Second, I don’t agree with the stackexchange answer you linked that the SAPT induction energy should be associated with the Coulomb potential. It’s not clear how to map the SAPT energies onto the nonbonded terms in a MM potential (i.e. Lennard-Jones plus Coulomb). The four groupings of energies reported at the end of a SAPT calculation (i.e. electrostatics, exchange, induction, and dispersion) are each composed of multiple terms, and the choice of how to assign these terms to the 4 groupings is not unique. For example, Psi4 reports the exchange-induction terms in the induction grouping and the exchange-dispersion terms in the dispersion grouping, but you could make an argument that it’s more appropriate to include these terms in the exchange grouping.

For these reasons, my recommendation is to fit all of the terms in the MM potential to the total SAPT interaction energy. The stackexchange answer that you linked clearly demonstrated that fitting only the Lennard-Jones term to the total SAPT energy gives a poor quality fit, but I think that fitting Lennard-Jones and Coulomb simultaneously to the total SAPT energy will be more accurate and more chemically rigorous than trying to fit only the Coulomb potential to the sum of the SAPT electrostatic and induction energies.

As a further complication, no matter what you decide to use as your fitting target, you may have a hard time getting accurate fit parameters if you’ve only varied the distance between the carbon atoms. You should also include some conformations that have the same carbon-carbon distance but different relative orientations between the methane groups. This will allow the fit to know which variations in energy come from the carbon-carbon interaction versus the other pairwise interactions between the molecules (i.e. carbon-hydrogen and hydrogen-hydrogen).

Thank you very much for your clarification @ccavender !

Before calculating the energies via SAPT, I optimized the geometry of the pairs by the MP2/aug-cc-pvtz method and obtained a set of bond lengths and angles relatively close to values ​​found in the literature.

I was trying to create a force field for the disyloxane molecule. Then, I obtained the Lennard-Jones parameters for H, O, C and Si through the energy curves of the hydrogen, water, methane and silane pairs, respectively. In order to obtain the Lennard-Jones parameters, I adjusted the SAPT energy curve of induction + dispersion as in the procedure performed on the first link I made available ( calculate-% 20lennard-jones-potential-with-quantum-mechanical-methods)

I found an article here on the forum ( commenting on the efficiency and performance of various levels of SAPT. Due to computational limitations, I choosed SAPT2+/aug-cc-pVDZ (second best method according to the article).

For the charges of the atoms of H, O, C, Si, I obtained them via ESP. The energy curves to obtain the intramolecular potential were obtained via HF/6-31g(2d,p) by varying the structure of the disyloxane molecule (this method was the one that best reproduced experimental data of bond lengths and angles)

Using the Lennard-Jones parameters for H, O, C, Si obtained by SAPT2+/aug-cc-pVDZ (induction + dispersion) in psi4, I did a molecular dynamics simulation in LAMMPS and obtained a density of ~ 0.78g/cm^3 vs 0.76 g/cm^3 experimental. I was “surprised” by the quality of this result.

However, I agree with what you said you should test more geometry settings for the pairs (I even saw in some articles the authors doing this). I think this will be the next step I will take!

I will try to apply your advice, if there are new results, I share.

Best regards,

For the charges of the atoms of H, O, C, Si, I obtained them via ESP.

This is fine. For the Lennard-Jones parameters, I would still recommend using the total SAPT energy as your fitting target instead of a subset of the SAPT energy groupings for the same reasons I gave above. Since you already have charges from ESP, you should fit the Lennard-Jones terms to the total SAPT energy minus the Coulomb energy from your ESP charges. This will allow you to use charges from ESP and LJ parameters from SAPT while ensuring that both nonbonded terms are consistent with the total SAPT energy.

I did a molecular dynamics simulation in LAMMPS and obtained a density of ~ 0.78g/cm^3 vs 0.76 g/cm^3 experimental.

This accuracy seems consistent with force fields for organic small molecules. If you think it is good enough for your application, then you can probably go with what you already have.