1-4 pairs and non bonded list

GROMACS version: 2019.6
GROMACS modification: No

Hi all,

I am currently running tests on reproducing the calculation of a protein in vacuum that was run using the amber14sb forcefield. I am writing the topology manually and comparing to the one where I run calling the forcefield itself. I know this seems redundant, but it will be important for my future work. I have a pair list defined with parameters and a non bonded section (nexc=3). Question is, are the 1-4 interactions counted twice with this set up? Or is it only taken from the [ pairs ] section since I have nexc=3? If the 1-4 is counted twice, then how do I make sure it is only taken into consideration once?

This is the topology file:
processedNa_f2sig05eps_explnb.top (145.9 KB)

Typically, force fields scale 1-4 interactions by some factor, so they have to be included, but not by conventional nonbonded interaction forms. So by setting nrexcl = 3, 1-4 interactions are removed (excluded). A [pairs] directive adds them back, but tells grompp to incorporate the scaling factors into the actual nonbonded form that mdrun subsequently uses. In this regard, a pair interaction is actually treated more like a bonded interaction than a nonbonded one.

Hi Justin,

thanks for the speedy reply. I think I am on the right track then. I have specified the non bonded interactions in [ nonbond_params] and all the 1-4 interactions separately in the [ pairs ] directive of my molecule, which I suppose will then add them back in, like you mentioned. So now the question is, because I need to have different fudge parameters, I have gen-pairs=no and I’m using this format in the pairs section, with V and W in sigma and epsilon form since func=2, right?

extra LJ or Coulomb [pairs] 2 2 fudge QQ (); qi; qj (e), V[7]; W[7]

I can see where to specify the different fudge QQ value, but how do I specify the different fudge LJ value? I’ve combed through the manual, but couldn’t find the answer to that one.


You don’t. That’s what [pairtypes] is for. Values there override those generated by the combination rule and default fudgeLJ (I think, please do a single-point energy calculation to confirm, I am not 100% sure if the default fudgeLJ gets applied or not).

Then wouldn’t that be the same as having everything explicitly in the [ pairs ] directive of the molecule as above but then giving epsilon values divided by two for fudge LJ=0.5? Given that V (r) = 4eps [(sig/r)^12-(sig/r)^6].

“The normal pair interactions, intended for 1-4 interactions, have function type 1. Function type 2 and the [ pairs_nb ] are intended for free-energy simulations” This is also written in the manual but in File formats — GROMACS 2021.4 documentation [ pairs ] have the possibility of having function type 1 or 2, whereas [ pairs_nb ] have function type 1. Is this correct?

Having never used those features, I can’t comment specifically, but type 2 and [pairs_nb] have specific use cases intended for applications related to free energy and corresponding topology transformations. I don’t know how they are interpreted in the context of standard MD simulations.

I think I have figured it out. Thank you for your replies, Justin. It’s been helpful to know what the program ignores and overrides when figuring out how to set up my topology.