Page 1 of 1

atom reordering with cutoffs

Posted: Wed Feb 03, 2016 3:15 pm
by brookscl
Peter et al.,

We are trying to finish up our CPHMD plugin to complement the recent implementation of the GBSW model in the CHARMM/OpenMM interface and are running into an issue that I think I ran across before but cannot find any detail of. This has to do with reordering of atoms to optimize the tile-based computations on the GPU (in cuda) when one has cutoffs. I am pretty sure this happens. We discovered an array (atomIndex) and a cuda platform kernel routine (CudaContext.cpp) where we believe this happens, our problem is 1) we want to gather the energies of interaction between specific sets of atoms and their neighbors at each timestep and 2) we use these energies to propagate the values of the charges on the atoms, and hence need to repopulate the charges with their updated values after the propagation step. It appears that this array, atomIndex, contains the information we need, but applying it to "reorder" our calculation doesn't seem to yield the correct result. I note that if we do no cutoff and don't use this array to reorder atoms, our results agree with those from or reference or CPU platform (i.e., CHARMM). Clearly we are not doing something correct, but we don't seem to be able to figure out where/why. Any help would be appreciated.

Thanks,

Charlie

Re: atom reordering with cutoffs

Posted: Wed Feb 03, 2016 3:46 pm
by peastman
Hi Charlie,

Two possible solutions, depending on exactly what you're doing.

First, you may simply be applying the atomIndex array in the wrong way. If atomIndex[ i] = j, that means the i'th atom in the internal ordering (as seen by GPU kernels) corresponds to atom j in the original ordering (as seen by public API methods). Whenever I get confused by this (which I do, often!), I look at the implementation of get/set routines for positions or velocities to make sure of the correct way to do the indexing: https://github.com/pandegroup/openmm/bl ... #L234-L291

The other possibility is simply to prevent any reordering that affects your force. We mostly try to do that, so we don't have to update lots of arrays of force field parameters after every reordering. The only type of reordering that can happen is to swap the indices of identical molecules. Mostly that means we reorder the waters. But since all the force field parameters are the same for every water molecule, that doesn't affect them.

The decision about what is "identical" is determined by the CudaForceInfo classes that all the forces register. They're described in more detail at http://docs.openmm.org/6.3.0/developerg ... -particles. So possibly you just need to define and register a CudaForceInfo. From your description, though, I'm not sure that's actually the right solution for your situation.

Peter