[gmx-developers] PBC with p2sub1 from pbc_dx() modifications

[David van der Spoel]

chris.neale at utoronto.ca wrote:
> Fantastic! Thank you so much for the details. It will take some weeks to test 
> that everything is working as it should, but I will be sure to post as soon as 
> everything is tested. In the end, I believe that a modern computer should still 
> get 0.5-1ns/day on a simple protein-membrane system, which is about 5 times 
> faster than charmm for me, making it a very worthwhile endevour. Once working, I 
> will be pulling a palmityl-(oleyl-phosphatidylethanolamine) out of a gram 
> negative outer menbrane protein (PagP) binding pocket by a derivative of replica 
> exchange and measuring the PMF making it quite important the the membrane is 
> happy. I really want to emphasize my appreciation for the rapid and in-depth 
> comments of Berk and Erik. And I won't hesitate to point out that the code is 
> generally far more readable than the humble Erik lets on.
> 

Just as a last comment, I'm not sure whether you mentioned this or not, 
but please base your efforts on the latest CVS, i.e. *not* the 3.3 tree. 
It will be much easier to merge with vanilla GROMACS later on, since 
Berk has already done quite a few modifications to the code there.

Good luck.

> Quoting Erik Lindahl <lindahl at sbc.su.se>:
> 
> 
>>Hi Chris,
>>
>>I've actually given this some thought in more general terms, and the  
>>appropriate long-term solution is probably to completely separate the  
>>asymmetric units/symmetry from the periodicity of the unit cell.
>>
>>
>>Let's start with the neighborsearching. This code might look  
>>complicated primarily due to a very deep nested loop over shifts and  
>>dimensions, but it's not really _that_ hard to understand. Referring  
>>to src/mdlib/ns.c:
>>
>>search_neighbors() is a wrapper function that does some setup,  
>>assigns grid cell indices to atoms (with fill_grid() in nsgrid.c) and  
>>then calls ns5_core() to do the real search work. All  
>>neighborsearching is based on charge group centers-of-geometry rather  
>>than atoms for reasons of efficiency.
>>
>>Still, _exclusions_ from NS need to be defined on an atom-by-atom  
>>basis, so we solve this in a pretty smart way with bit masks.  Charge  
>>groups are currently allowed to contain up to 32 atoms, so before  
>>doing advanced grid stuff and locating neighbors for a particular  
>>group we go through the (usually very short) list of excluded atoms  
>>for this ("i particule") group and for each of these ("j particle")  
>>set the bits of the corresponding excluded atoms to 1. You probably  
>>won't need to alter this, it's just to give you a hint what it does.
>>
>>in ns5_core() we go through the neighbor grid cells that might  
>>contain neighbor particles - the number of cells we need to check  
>>depends on the box geometry, so we first determine exactly how many  
>>neighbor cells we should check, and in which directions. However,  
>>instead of checking periodic boundary conditions for every single  
>>distance we use "shift vectors" here too. Think of an imaginary 1D- 
>>box where the current cell is the leftmost, and we want to find  
>>neighbors both to the right (next box), and to the left (using PBC,  
>>in the rightmost box). This corresponds to a list of two shift  
>>vectors. First we search the next box without doing anything special.  
>>For the second step, we first add a shift vector corresponding to the  
>>1D-box length, which you can think of as moving the entire reference  
>>box to the right of the system, and then search the rightmost box  
>>without any need for checking PBC.
>>
>>You WILL need to understand this "nonbonded shift" concept - there  
>>are some illustrations in the manual that might help, as well as the  
>>Gromacs papers!
>>
>>Finally, we check distances between charge groups and eventually go  
>>through atoms and put them in the lists. This is just a bunch of  
>>special cases to optimize the neighborlists as much as possible, e.g.  
>>VdW-only, or coulomb-interactions-only between pairs of water molecules.
>>
>>As explained in the manual, we need to store the shift indices/ 
>>vectors for each neighborlist to be able to calculate the virial as a  
>>single sum after calling the nonbonded forces instead of doing it  
>>inside the kernel. I somewhat suspect that it might be easiest for  
>>you just to copy coordinates to get a full periodic cell so you don't  
>>have to worry about rotations in the shift - that could namely affect  
>>the innerloops too, and you don't want to hack assembly ;-)
>>
>>
>>Now for something completely different: shifts.
>>
>>
>>There is an related, but slightly different, periodicity problem with  
>>bonded interactions. We really should have used a different name, but  
>>unfortunately it's called "molecular shift" or even only "shift" in  
>>the code. So sorry. Mea culpa...
>>
>>The graph code goes through all molecules at setup time and makes  
>>them continuous, i.e. removes PBC. We then re-apply PBC, but keep  
>>track of how each atom has been moved/shifted in terms of integer  
>>multiples of the unit cell vectors.
>>
>>Neighborsearching and nonbonded force calculation normally puts all  
>>charge groups inside the central unit cell, but the molecular shifts  
>>enable us to rapidly make all molecules continuous again without  
>>evaluating a single conditional. I'm not sure if you need bother with  
>>this, but I think we explained this too in more detail in the latest  
>>JCC Gromacs paper.
>>
>>PME is fairly independent from everything else, in src/mdlib/pme.c.  
>>You provide coordinates and triclinic box vectors and get potential,  
>>forces on each atom, and long-range virial back.
>>
>>My recommendation would be to start with neighborsearch for a trivial  
>>system without bonded interactions or exclusions, and make sure you  
>>get the neighborlist correct :-)
>>
>>Cheers,
>>
>>Erik
>>
> 
> _______________________________________________
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the 
> www interface or send it to gmx-developers-request at gromacs.org.


-- 
David.
________________________________________________________________________
David van der Spoel, PhD, Assoc. Prof., Molecular Biophysics group,
Dept. of Cell and Molecular Biology, Uppsala University.
Husargatan 3, Box 596,  	75124 Uppsala, Sweden
phone:	46 18 471 4205		fax: 46 18 511 755
spoel at xray.bmc.uu.se	spoel at gromacs.org   http://folding.bmc.uu.se
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++