[gmx-users] Transition difficulties: version 3.3.3 to 4.0.3 regarding pull_geometry=distance

Steve Fiedler fiedler at umich.edu
Wed Jan 28 06:52:06 CET 2009 

Hi,

The permeants in my bilayer system are improperly pulled in calculations 
using all versions of Gromacs 4.0.x including the VERSION 
4.0.99_development_20090120.  Since, due to restrictions, I can not 
provide the topology file for my system to bugzilla, I constructed a 
similar prototypical system using Prof. Tieleman's POPC bilayer 
coordinates and field parameters.  A buckyball permenat doped in the 
center of this system did indeed react similarly in a 4.0.3 calculation 
with an improper pull away from the equilibrated 3.3.1 constraint force 
coordinates.  It was heartening to see that this problem was fixed for 
that system, by using VERSION 4.0.99_development_20090120.  
Unfortunately, again the problem remained for my actual system, even 
with use of the newest version of code.

Since the cylinder option is now operational, I have a viable 
work-around for this problem.  I am quite appreciative of the attention 
that this issue has received.  If there is continued interest, and a 
willingness for an engineer to receive the topology file in a less 
public forum than bugzilla, I believe we could proceed with the 
debugging process.

Thank you,

Steve Fiedler




Berk Hess wrote:
> Hi,
>
> Just to be sure, the pull code is producing the correct results now?
> The only problems should be in several analysis tools.
>
> The analysis tools have been written only with analysis of molecules
> or parts of molecules in mind. When you want to analyze groups
> that cover two molecules or more there will be pbc problems.
> These are not simple to fix. If the distances between all the atoms
> in the group are less than half a box length there is a unique COM.
> The procedure should be similar to what trjconv -pbc cluster does.
>
> I can try to implement this for g_traj and g_dist.
> Are there other tools that cause problems?
>
> Berk
>
> > Date: Mon, 26 Jan 2009 17:46:26 -0500
> > From: fiedler at umich.edu
> > To: gmx-users at gromacs.org
> > Subject: Re: [gmx-users] Transition difficulties: version 3.3.3 to 
> 4.0.3 regarding pull_geometry=distance
> >
> > Hi Chris,
> >
> > Your suggestions were helpful. Due to the nature of this problem, there
> > were limitations in the application of the sequence of steps that you
> > recommended. Specifically, the difference in the equilibrated structure
> > from version 3.3.x to the calculated reference center of mass from
> > version 4.0.x causes an unreasonable "jump" with constraint force
> > calculation. Resultant bad contacts prohibit generation of a pull
> > output file.
> >
> > More generally though, I believe I understand your broader point: the
> > sensitive nature of the periodic box with respect to the pull code and
> > possible discrepancies in the pbc treatment by various utilities, may
> > cause confusion with the diagnostic process. The test bilayer system I
> > created however, centered distantly from the x and y box edges, may
> > effectively remove this set of concerns. For example, this
> > configuration is relatively unchanged upon "recentering", however it is
> > still susceptible to the problem that is the subject of this discussion.
> >
> > If I can reproduce this phenomenon with a publicly available topology
> > file, I can provide the necessary input files for inspection.
> >
> > Thank you again,
> >
> > Steve Fiedler
> >
> >
> >
> > chris.neale at utoronto.ca wrote:
> > > Hi Steve,
> > >
> > > what I intended to suggest was actually something different (and much
> > > easier).
> > >
> > > The idea is not that you need some special system to be able to
> > > utilize the pull code, but that the pull code is correct whereas the
> > > g_dist and g_traj programs are not as good at treating pbc in the way
> > > that one desires.
> > >
> > > I suggest the following.
> > >
> > > 1. Take your original system and run the pull code for a very short
> > > simulation. Use the last line of the output to calculate the relevant
> > > displacement
> > >
> > > 2. Now use trjconv -b -e to get the last frame of the .xtc that
> > > resulted from that short MD run as a .gro file, call it final.gro. I
> > > suspect that your groups are not entirely in the same simulation box
> > > in final.gro.
> > >
> > > 3. Now make a new .ndx file from that .gro and give it a single
> > > residue that is near your binding pocket, call it R_1
> > >
> > > 4. Now apply trjconv -center -pbc mol -ur compact while selecting R_1
> > > for centering, call the new .gro file final_center.gro
> > >
> > > 5. Visualize final_center.gro and ensure that all of your relevant
> > > atoms are in the same image in the way that puts the minimum distance
> > > between them along a path that is entirely contained within the unit
> > > cell. If not, go back to step 3 and try making a group R_2, etc.
> > > until this process works. NOTE: you might think that giving trjconv
> > > -center the relevant groups that you use for pulling will be a good
> > > idea here, but it is not. The problem there is that the atoms may be
> > > "centered" by placing half on the left boundary and half on the right
> > > boundary. I find using one logically selected residue or atom is the
> > > best method here.
> > >
> > > 6. Assuming that you got what you wanted in step 5, now run g_traj 
> and
> > > g_dist on final_center.gro. In my case, I found that g_traj and 
> g_dist
> > > give the same answer as the pull code output when I am using
> > > final_center.gro, but not always when I am using final.gro.
> > >
> > > *** I always laugh when these problems arise because, in an important
> > > sense, the protein *did* jump out of the simulation box... at 
> least as
> > > far as g_traj and g_dist are concerned. This, we must hope, is
> > > correctly treated in the pull code even though it is incorrectly (or
> > > at least unintuitively) treated by g_traj and g_dist.
> > >
> > > Chris.
> > >
> > > -- original message --
> > >
> > > Hi,
> > >
> > > Thank you Berk and Chris for the suggestions.
> > >
> > > To address the possibility that this issue is related to periodic
> > > boundaries, I used two approaches:
> > > 1. The pull group of interest (permeant) was centered in the x-y plane
> > > of the box using Chris' approach. I then used the genconf utility to
> > > replicate my lipid box to a 9x9 grid in the x-y plane and removed all
> > > but the center box. This generated the coordinates for a bilayer 
> system
> > > with all lipid molecules inside a box and intact. The discrepancy
> > > between the grompp (version 4.0.3) output and distances as 
> calculated by
> > > g_traj (version 4.0.3) persist, 2.667 vs. 0.3996 nm.
> > > 2. I constructed a three atom system containing 2 reference atoms of
> > > type A, and a "pull" atom of type B. Proper output from grompp was
> > > observed for all coordinates of both the reference and pulled atoms,
> > > include coordinates for atoms moved outside the box in the x-y plane.
> > > The coordinate, topology, and run control parameter file are given 
> below.
> > >
> > > If there are additional suggestions, I would be greatly appreciative.
> > >
> > > Thank you,
> > >
> > > Steve Fiedler
> > >
> > > -----------------
> > > conf.gro
> > > Three atoms
> > > 3
> > > 1AAA A 1 1.500 1.500 1.000
> > > 2AAA A 2 0.500 1.500 1.000
> > > 3BBB B 3 -1.500 1.500 1.700
> > > 3.00000 3.00000 3.00000
> > > -----------------
> > > index.ndx
> > > [ System ]
> > > 1 2 3
> > > [ Ref ]
> > > 1 2
> > > [ Pulled ]
> > > 3
> > > -----------------
> > > grompp.mdp
> > > title = ThreeAtoms
> > > integrator = md
> > > dt = 0.001
> > > nsteps = 1
> > > ns_type = grid
> > > pbc = xyz
> > > coulombtype = shift
> > > rlist = 1.4
> > > rcoulomb = 1.4
> > > rvdw = 1.4
> > > tcoupl = no
> > > pcoupl = no
> > > constraint_algorithm = shake
> > > shake_tol = 1e-4
> > > gen-vel = no
> > > gen-temp = 0
> > >
> > > nstxout = 1
> > > nstvout = 0
> > > nstfout = 0
> > >
> > > pull = umbrella
> > > pull_geometry = distance
> > > pull_dim = N N Y
> > > pull_start = no
> > > pull_init1 = 0.7
> > > pull_group0 = Ref
> > > pull_group1 = Pulled
> > > pull_k1 = 10000
> > > -----------------
> > > topology.top
> > > ; topology for two partially charged atoms
> > >
> > > [ defaults ]
> > > ; nbfunc comb-rule gen-pairs fudgeLJ fudgeQQ
> > > 1 3 yes 0.125 0.5
> > >
> > > [ atomtypes ]
> > > ;name mass charge ptype sig eps
> > > A 1000.0000 0.000 A 0.50000 9.90000
> > > B 9.0000 0.000 A 0.30000 9.00000
> > >
> > > [ nonbond_params ]
> > > ; i j func sig eps
> > >
> > > [ moleculetype ]
> > > AAAA 1
> > >
> > > [ atoms ]
> > > ; nr type resnr residue atom cgnr charge mass
> > > 1 A 1 AAA A 1 0.000 1000.0000
> > >
> > > [ moleculetype ]
> > > BBBB 1
> > >
> > > [ atoms ]
> > > ; nr type resnr residue atom cgnr charge mass
> > > 1 B 1 BBB B 1 0.000 9.00
> > >
> > > [ system ]
> > > ; name
> > > Three atoms
> > >
> > > [ molecules ]
> > > ; name number
> > > AAAA 2
> > > BBBB 1
> > >
> > >
> > >
> > >
> > >
> > > Chris Neale wrote:
> > >> I just checked similar simulations of mine and Berk's suggestion
> > >> accounts for similar discrepancies that I notice on a quick
> > >> evaluation where g_traj and g_dist fail to give me the same distance
> > >> as I obtain from the pull pos.xvg file. As Berk suggests, once I
> > >> first trjconv -center -pbc mol -ur compact (giving an appropriate
> > >> residue for centering that puts all relevant pulled atoms in the 
> same
> > >> box) then g_traj and g_dist both give me the exact same answer as I
> > >> calculate based on pull pos.xvg. Chris -- original message -- Hi,
> > >> There could be a problem with periodic boundary conditions. Do you
> > >> have multiple molecules in a pull group, or broken molecules? In 
> that
> > >> case the COM position of 3.3.3 and g_traj are both incorrect. The
> > >> pull code in 4.0 grompp and mdrun are (as far as I know) always
> > >> correct. Berk
> > >
> > > _______________________________________________
> > > gmx-users mailing list gmx-users at gromacs.org
> > > http://www.gromacs.org/mailman/listinfo/gmx-users
> > > Please search the archive at http://www.gromacs.org/search before
> > > posting!
> > > Please don't post (un)subscribe requests to the list. Use thewww
> > > interface or send it to gmx-users-request at gromacs.org.
> > > Can't post? Read http://www.gromacs.org/mailing_lists/users.php
> > >
> > >
> >
> >
> > --
> > Steve Fiedler, Ph.D.
> > Research Fellow
> > Department of Mechanical Engineering
> > The University of Michigan
> > 2024 G.G. Brown
> > 2350 Hayward St.
> > Ann Arbor, MI 48109-2125
> >
> >
> > _______________________________________________
> > gmx-users mailing list gmx-users at gromacs.org
> > http://www.gromacs.org/mailman/listinfo/gmx-users
> > Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
>
> ------------------------------------------------------------------------
> What can you do with the new Windows Live? Find out 
> <http://www.microsoft.com/windows/windowslive/default.aspx>
> ------------------------------------------------------------------------
>
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php

More information about the gromacs.org_gmx-users mailing list