[gmx-developers] Extracting i-j forces from nb_generic.c

James jamesresearching at gmail.com
Thu Jan 9 07:43:26 CET 2014

Dear Erik, David and other Gromacs developers,

Thank you very much for your replies. They were very helpful! I have run
with just one thread and I can now see how to extract i-j VDW and
short-range coulomb potentials (shared equally between pairs) and by
tracing the temperature evaluation I think I have found the appropriate
place for extracting kinetic energy and velocities.

I wonder if I may follow up about the atomic numbering. I will need to move
to parallel calculation as soon as possible, and I have been thinking about
the comments:

(Erik) "you can write functionality to move the properties you need (just
as we move charge, etc) into the local topology when we do rebalancing
between domains"

I'm not sure what exactly is meant by "local topology". Where is
"rebalancing between domains" performed? Would you mind expanding this
explanation a little more or giving an example?

(David) "There is a function in mtop_util.h that will give you the global
atom number"

Are you referring to the gmx_mtop_atomlookup_init function? I looked at
other places in the code where this is called; I assume I would have to
somehow pass *mtop to the force subroutine? If evaluating this is
expensive, as Erik suggests, I suppose I could do this outside of the force

Thank you again for your help.

With best regards,


On 31 December 2013 18:15, Erik Lindahl <erik.lindahl at scilifelab.se> wrote:

> Hi James,
> On 31 Dec 2013, at 08:15, James <jamesresearching at gmail.com> wrote:
> > Dear Gromacs developers,
> >
> > To measure thermal conductivity using Green-Kubo I am trying to extract
> the i-j atomic forces from nb_generic.c and I wonder if anyone could help
> me with some basic questions?
> >
> > 1. If I understand correctly, ii and jj are atom numbers for atom i and
> atom j. These don't seem to correspond to the atoms in the order that they
> were read in (PDB file), however. How can I recover this order, or what is
> the new order?
> That would be “ii” and “jnr” in the master branch. If you are only using a
> single thread this should correspond to the atom order in the input
> (although C always starts counting at 0). However, when we run in parallel
> the domain decomposition will mean that the particle indices change
> (frequently) as particules move across node boundaries.
> There are routines you can use to find the global atom indices, but you
> don’t want to call those from the innermost kernel. Instead I would
> recommend using a single thread while you develop, and later when you want
> to enable parallelization you can write functionality to move the
> properties you need (just as we move charge, etc) into the local topology
> when we do rebalancing between domains.
> > 2. What units are the positions ix, iy, iz (again in nb_generic.c)?
> nm - same as in the gro files, but it differs from PDB files that use
> Angstrom.
> >
> > 3. nb_generic is within the loop i = i0 to i1-1 in nonbonded.c. What is
> this loop? Somehow related to different interaction types?
> At one point we experimented with multithreading parallelism over this
> loop (thus the indices), but now the outermost loop is over different
> neighborlists. A certain particle will have different neighborlists for
> neighbors that only have charge, only Lennard-Jones, and both LJ+charge.
> Having said that, for the future we are likely moving entirely to the new
> style “verlet” kernels currently used for GPUs (and that will likely be
> default for CPUs in Gromacs-5.0 too), but we still need to write a generic
> kernel there.
> > 4. Later, I will also want to extract the instantaneous kinetic and
> potential energy of the atom in the same time-step - any advice is very
> welcome.
> For kinetic energy you can just look at how we calculate temperature -
> this is easy. Potential energy is complex if you use lattice summation
> (PME), since you then would need to modify the grid interpolation to get
> potentials of individual particles.
> >
> > 5. Out of interest, I notice many of the files in the release-4-6 branch
> are in include/ but the same files are in legacyheaders/ in the master
> branch. Why is this? Which should I be developing against?
> We are moving pretty fast right now due to a switch from C to partial C++,
> which also involves a lot of cleaning-up and modularization. The master
> branch is the future, but legacyheaders is stuff that still hasn’t been
> modularized.
> Cheers,
> Erik
> --
> Gromacs Developers mailing list
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List before
> posting!
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developersor send a mail to
> gmx-developers-request at gromacs.org.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-developers/attachments/20140109/b035381f/attachment.html>

More information about the gromacs.org_gmx-developers mailing list