[gmx-users] Re: Implementing Radial Distance Constraints

Mark Abraham Mark.Abraham at anu.edu.au
Tue Apr 18 03:39:43 CEST 2006

Bob Johnson wrote:
>>please be more clear.
>>if we assume your cnt is centered along the z-axis, given roughly by
>>x^2+y^2 = d0. Am I correct in assuming that you want all other particles
>>to have x^2+y^2 <= d1, where d1 < d0 ?
> The molecules in my system adsorb to the tube's outer wall. They are not inside 
> the tube. Typically, the adsorbates reside 3.5A above the nanotubes surface (a
> typical distance for molecules interacting via VdW with an aromatic molecule).
> I want to constrain molecules to always adopt a monolayer about the nanotube.
> Therefore, I want to put a constraint on certain atoms to reside 3.5A above the
> nanotube. Let's say my nanotube has a radius of 0.5A. I would then want to
> implement a potential given by V = 1/2 k (R-0.5A)^2 where R = sqrt(x^2+y^2).
> This way, the atoms are able to slide along the z and cirumferential
> directions.

That potential will attempt to drag the adsorbates down to the radius of 
the nanotube - vdW forces will oppose this and equilibrium will occur 
somewhere, but probably inside 3.5A. I think there are less perturbative 
ways of applying these constraints. A potential like half of Figure 
4.1.3 in the gromacs manual would make more sense, only switching on 
when R > 3.5+0.5+x for some small x that acknowledges that 3.5 is an 
average distance. That is

V = 0 for R <= R0
V = 0.5 k ( sqrt(x^2+y^2) - R0 ) ^2 for R > R0

Since you expect R-R0 to be small, you might consider implementing the 
approximate potential V = 0.5 k ( R^2 - R0^2 ) for computational efficiency.

As far as finding out what to change goes, you'll need to change grompp 
to recognize a new type of cylindrical constraint, and then provide a 
function that evaluates the potential. There are a bunch of directives 
that affect the function of distance and orientation restraints for MD 
which won't want to concern you, so I suggest using the code for angle 
restraints as a template to adapt to your purposes. In the gromacs src 
directory, "grep angle_restraints */*.[ch]" and also with "*/*/*.[ch]" 
to find the headers and code sections where grompp interprets the .top 
file. That will give a clue for the names of relevant data structures 
and what needs to be added where. You will need to modify 
src/gmxlib/bondfree.c to add a new potential type - grepping for 
"angres" will help you find all the places you need to add the new 
function and the supporting structure to work out when it gets called.


More information about the gromacs.org_gmx-users mailing list