[gmx-users] Re: What is the purpose of the [ pairs ] section?

[Bogdan Costescu]

Hi Andrew!

I've also been a bit confused by the names when I started looking at
the force field definitions in GROMACS. Although Justin's replies
contain all the information (as usual :)), let's try to shed some
light with an example: consider a chain formed by atoms i-j-k-l-m-n; -
means that there is a bond between atoms.

nrexcl specifies how many neighbors should be excluded when building
lists of nonbonded interactions. nrexcl = 3 means that the first,
second and third order neighbors (j, k and l for atom i; k, l and m
for atom j; l, m and n for atom k; m and n for atom l; n for atom m)
are excluded, such that the following pairs will be put on the lists
to be passed to the nonbonded kernels:

i m
i n
j n

For another force field with nrexcl = 2, the following pairs would be
passed to the nonbonded kernels:

i l
i m
i n
j m
j n
k n

For another force field with nrexcl = 5, no pairs would be passed to
the nonbonded kernels. More generally, any pairs which are neighbors
of higher order than nrexcl will be put on the nonbonded lists.
Parameters for these nonbonded interactions will be automatically
computed from the topology data, based on the combination rules. If
you have empty [bonds], [angles], [dihedrals] and [pairs] sections in
the topology file, only the above interactions are considered for
potential energy calculations. These nonbonded lists are dynamical and
are remade based on neighbor search frequency; when using cutoffs, a
pair whose atoms are further apart than the cutoff is not added to the
list.

Bonded interactions are specified in the [bonds], [angles] and
[dihedrals] sections. Parameters for the bonded interactions are taken
directly from the topology data after matching atom types. Any of
these sections can be empty, if the force field is so designed; for
OPLS-AA all these sections are present. There can be overlaps, f.e. if
nrexcl = 2 and the topology file contains:

[dihedrals]
i j k l
j k l m
k l m n

then atoms i and l, j and m and k and n will interact both as part of
the dihedral angle potential and as nonbonded pair. Some force fields
(like OPLS-AA) specify that the 1-4 interactions ((i l), (j m) and (k
n) in the above example) should indeed be treated both as part of a
dihedral angle potential and as nonbonded pairs; when treating them as
nonbonded pairs there might be a scaling factor (0.5 for OPLS-AA).
When there is no scaling factor (in other words, scaling factor is 1),
a definition like above with nrexcl = 2 is all that is needed. For a
scaling factor different from 1, the GROMACS implementation uses an
explicit pairs list, called [pairs]. All 1-4 interactions should be
included in the [pairs] section, which for the above example would
look like:

[pairs]
i l
j m
k n

It's the user's responsibility to add all needed (1-4 in case of
OPLS-AA) pairs to the [pairs] section; pdb2gmx normally does this, but
if you build the topology through other means then you have to make
sure these pairs are generated. The interactions for these explicit
pairs are always computed - even though they are computed through the
nonbonded kernels, they are not subject to cutoffs. These explicit
pairs also need parameters. The parameters are automatically generated
when gen-pairs = yes from the topology data, based on the combination
rules and the fudge (=scaling) factors. If gen-pairs = no, then the
[pairs] section should contain both the atom indeces and the
parameters. I found the name gen-pairs a bit misleading, I would have
called it gen-pairs-params or something similar as it is only taken
into consideration when looking for parameters; the [pairs] list is
not automatically updated when generating the .tpr file. The [pairs]
section offers some flexibility, as not only the 1-4 interactions but
any other interaction that should be treated the same way can be
included; a possible limitation is that only one fudge factor can be
set for the whole list of explicit pairs.

Care should be taken to not calculate the same interaction twice. With
nrexcl = 3 and:

[pairs]
i m

the nonbonded interaction between atoms i and m will be computed
twice: once as part of explicit pairs list and once as part of
nonbonded lists.

Finally, I would like to point out that most of these can be found by
reading the manual and experimenting with a toy molecular system (like
the example above) using gmxdump and a combination of mdrun -rerun and
g_energy.

Cheers,
Bogdan

On Mon, Dec 10, 2012 at 9:13 PM, Andrew DeYoung <adeyoung at andrew.cmu.edu> wrote:
> Javier and Mark,
>
> Thank you SO much!  That is so helpful to me.
>
> If you have time, I have two follow-up questions:
>
> (1) From your answer, it seems that in OPLS, 1-4 interactions will be taken
> into account (scaled by fudge) ONLY if the 1-4 atomic indices are specified
> in [ pairs ], even if gen-pairs = yes and sigma and epsilon values are
> provided in [ atomtypes ] in ffnonbonded.itp.
>
> But what about 1-5, 1-6, 1-7, ..., 1-Infinity nonbonded interactions, in a
> protein for example?  I think that in OPLS, such interactions are to be
> considered at full strength (no scaling down).  My question is: are these
> 1-5, 1-6, 1-7, ..., 1-Infinity nonbonded interactions automatically taken
> into account if gen-pairs = yes?
>
> (2) Suppose that I run grompp (where 1-4 interactions are specified in [
> pairs ] and gen-pairs = yes) and grompp gives me the following output:
>
> Generated 276 of the 276 non-bonded parameter combinations
> Generating 1-4 interactions: fudge = 0.5
> Generated 276 of the 276 1-4 parameter combinations
>
> Does this output indeed imply that the ONLY nonbonded interactions
> considered are 1-4 interactions?  In other words, does this output suggest
> that NO 1-5, 1-6, 1-7, ..., 1-Infinity nonbonded interactions are
> considered?  If so, why is this the case if I have set gen-pairs = yes?
>
> Andrew DeYoung
> Carnegie Mellon University
>
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> * Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists