Energy groups find self interacting Na+ and Cl+ for Coulombic interaction energy?

User discussions
1

GROMACS version:2022v2
GROMACS modification: No
Hi,
I am trying to calculate the Coulombic interaction energy between a single Na+ and Cl- ions. Using experimental distance of 0.276 nm, I can calculate 503 kJ/mol from the simple coulombic potential energy equation.

When I try to calculate the same by Gromacs using energy groups, it gives me the following output.
Energy Average Err.Est. RMSD Tot-Drift

Coulomb (SR) -563.928 – 0 0 (kJ/mol)
Coul. recip. 60.4996 – 0 0 (kJ/mol)
Potential -489.701 – 0 0 (kJ/mol)
Coul-SR:S1-S1 -204.026 – 0 0 (kJ/mol)
LJ-SR:S1-S1 0 – 0 0 (kJ/mol)
Coul-SR:S1-S2 -155.875 – 0 0 (kJ/mol)
LJ-SR:S1-S2 13.7274 – 0 0 (kJ/mol)
Coul-SR:S2-S2 -204.026 – 0 0 (kJ/mol)
LJ-SR:S2-S2 0 – 0 0 (kJ/mol)

So it looks like when creating energy groups, it finds self interaction terms of Na+Na+ and Cl-
I am using Verlet cutoff scheme (unfortunately group option is deprecated). Is it interacting with periodic images? Since the new versions do not allow other than Verlet cutoff scheme and Verlet only works with pbc, I could not test it on 2022 v 2.

But I have installed Gromacs 2018 and tried on that with pbc=no and cutoff-scheme=group with all cutoffs and nstlist stting to 0. It does not produce any artificial interaction and produces the interaction correctly:
Energy Average Err.Est. RMSD Tot-Drift

LJ (SR) 13.7274 – 0 0 (kJ/mol)
Coulomb (SR) -503.389 – 0 0 (kJ/mol)
Potential -489.662 – 0 0 (kJ/mol)
Coul-SR:S1-S1 0 – 0 0 (kJ/mol)
LJ-SR:S1-S1 0 – 0 0 (kJ/mol)
Coul-SR:S1-S2 -503.389 – 0 0 (kJ/mol)
LJ-SR:S1-S2 13.7274 – 0 0 (kJ/mol)
Coul-SR:S2-S2 0 – 0 0 (kJ/mol)
LJ-SR:S2-S2 0 – 0 0 (kJ/mol)

2

You seem to be using PME in the first case and not in the second, so these are completely different calculations. PME computes periodic Coulomb interactions.

3

Thanks for your reply. That is my bad. Yes, you are right, I used pme. Using cutoff instead mostly fixed the problem. But, i still get artificial self-self interactions unless I give really long cutoffs. I tried 180 nm but still small self interaction.
Energy Average Err.Est. RMSD Tot-Drift ------------------------------------------------------------------------------- Coulomb (SR) -503.391 – 0 0 (kJ/mol) Potential -489.663 – 0 0 (kJ/mol) Coul-SR:S1-S1 -3.838 – 0 0 (kJ/mol) LJ-SR:S1-S1 0 – 0 0 (kJ/mol) Coul-SR:S1-S2 -495.715 – 0 0 (kJ/mol) LJ-SR:S1-S2 13.728 – 0 0 (kJ/mol) Coul-SR:S2-S2 -3.838 – 0 0 (kJ/mol) LJ-SR:S2-S2 0 – 0 0 (kJ/mol)

4

I think you should get what you want when you set coulomb-modifier=none

5

Thanks a lot. That solved my problem.
I am currently working on charged vs non charged groups and would like to get electrostatic interaction of groups as good as possible. So using pme, and any switching is not good for at least being on the same page with my own script, Amber or other software comparsion. I see you had already mention that point on the mdp options help page :“coulomb-modifier=none: This can be useful when comparing energies with those computed with other software.”

Thanks a lot for your help.