[gmx-users] restraining atoms to the plane at the bilayer center: pull code ?
maria goranovic
mariagoranovic at gmail.com
Thu Sep 9 10:54:01 CEST 2010
Hmm. I am getting the desired output now, at least as far as the pulled
atoms are concerned, they are all close to the bilayer center (was not so in
a free simulation), and therefore my suspicion is that the code is working.
Just to re-iterate, the atoms to be pulled to the bilayer center were
present at the lipid-water interface in both leaflets.
Here is the (apparently successful) pull code:
; COM PULLING
pull = umbrella
pull_geometry = position
pull_dim = N N Y
pull_start = no
pull_nstxout = 5000
pull_nstfout = 5000
pull_ngroups = 2
pull_group0 = POPC
; global atom number of a POPC atom in the middle of the xy box
pull_pbcatom0 = 3100
pull_group1 = C35TOP
; global atom number of a C35 atom in the middle of the xy box
pull_pbcatom1 = 5467
pull_init1 = 0 0 0.0
pull_rate1 = 0
pull_k1 = 500
pull_vec1 = 0 0 0
pull_group2 = C35BOT
; global atom number of a C35 atom in the middle of the xy box
pull_pbcatom2 = 7755
pull_init2 = 0 0 0.0
pull_rate2 = 0
pull_k2 = 500
pull_vec2 = 0 0 0
##########################
Here is the output from grompp
Pull group natoms pbc atom distance at start reference at t=0
0 4992 3100
1 16 5467 0.000 0.000 3.184 0.000 0.000 0.000
2 16 7755 -0.000 0.000 -3.160 0.000 0.000 0.000
###########################
The average COM pull En. and potential from the log file after 100 ps.
2.96054e+03
###########################
And here is the output in the log files, detailing the pulling parameters
pull = umbrella
pull_geometry = position
pull_dim (3):
pull_dim[0]=0
pull_dim[1]=0
pull_dim[2]=1
pull_r1 = 1
pull_r0 = 1.5
pull_constr_tol = 1e-06
pull_nstxout = 5000
pull_nstfout = 5000
pull_ngrp = 2
pull_group 0:
atom (4992):
atom[0,...,4991] = {0,...,4991}
weight: not available
pbcatom = 3099
vec (3):
vec[0]= 0.00000e+00
vec[1]= 0.00000e+00
vec[2]= 0.00000e+00
init (3):
init[0]= 0.00000e+00
init[1]= 0.00000e+00
init[2]= 0.00000e+00
rate = 0
k = 0
kB = 0
pull_group 1:
atom (16):
atom[0]=5026
atom[1]=5114
atom[2]=5202
atom[3]=5290
atom[4]=5378
atom[5]=5466
atom[6]=5554
atom[7]=5642
atom[8]=5730
atom[9]=5818
atom[10]=5906
atom[11]=5994
atom[12]=6082
atom[13]=6169
atom[14]=6258
atom[15]=6346
weight: not available
pbcatom = 5466
vec (3):
vec[0]= 0.00000e+00
vec[1]= 0.00000e+00
vec[2]= 0.00000e+00
init (3):
init[0]= 0.00000e+00
init[1]= 0.00000e+00
init[2]= 0.00000e+00
rate = 0
k = 500
kB = 500
pull_group 2:
atom (16):
atom[0]=6434
atom[1]=6522
atom[2]=6610
atom[3]=6698
atom[4]=6786
atom[5]=6874
atom[6]=6962
atom[7]=7050
atom[8]=7138
atom[9]=7226
atom[10]=7314
atom[11]=7402
atom[12]=7490
atom[13]=7578
atom[14]=7666
atom[15]=7754
weight: not available
pbcatom = 7754
vec (3):
vec[0]= 0.00000e+00
vec[1]= 0.00000e+00
vec[2]= 0.00000e+00
init (3):
init[0]= 0.00000e+00
init[1]= 0.00000e+00
init[2]= 0.00000e+00
rate = 0
k = 500
kB = 500
#########################
I would want to make sure I am doing this right
Best Wishes
Maria
On Wed, Sep 8, 2010 at 4:09 PM, <chris.neale at utoronto.ca> wrote:
> Maria,
>
> I highly doubt that you are correct. If you want more help, please copy and
> paste your commands and some relevant output back to the list.
>
> Chris.
>
> --original message --
>
> Message: 4
> Date: Wed, 8 Sep 2010 15:07:15 +0200
> From: maria goranovic <mariagoranovic at gmail.com>
> Subject: Re: [gmx-users] restraining atoms to the plane at the bilayer
> center: pull code ?
> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
> Message-ID:
> <aANLkTi=RQv+Ozsff0u5pMXsUAq-jpJaGhnNVpDrJM=ZY at mail.gmail.com>
> Content-Type: text/plain; charset="iso-8859-1"
>
> Dear Chris,
>
> Thank you for the code. I did check out pull_pbcatomN, and thank you for
> the
> heads up that it is a global index. I used a central atom in the box for
> each of the two groups as pull_pbcatom.
>
> After using the code below with a rate constant of ~ 5500 for about 100 ps,
> I realized that the atoms in different leaflets were being pulled in the
> same direction. I wanted certain atoms of either leaflet to be restrained
> to the bilayer center. However, the code pulls all the atoms in the same
> direction, such that atoms from one leaflet do approach the bilayer center,
> but those from the other leaflet tend to get further away from the bilayer
> center. Thus, one actually needs two groups to be pulled, one each for the
> atoms in each leaflet. Once this was fixed, the code worked just fine.
>
> Again, thank you.
>
> Maria
>
>
>
> On Tue, Sep 7, 2010 at 3:40 PM, <chris.neale at utoronto.ca> wrote:
>
> [Hide Quoted Text]
> Maria, try this. There actually is a lot of this on the mailing list, so I
> suggest checking it a little deeper for your next querry, or at least
> outlining how you looked and what you found so that it is clear you have
> tried.
>
> Also, read about pull_pbcatomN and think carefully about how you want to
> set that up. It is a *global* index.
>
> ; COM PULLING
> pull = umbrella
> pull_geometry = position
> pull_dim = N N Y
> pull_start = no
> pull_nstxout = 500
> pull_nstfout = 500
> pull_ngroups = 1
> pull_group0 = <<<bilayer-group>>>
> pull_pbcatom0 = 0
> pull_group1 = <<<pulled-group>>>
> pull_pbcatom1 = 0
> pull_init1 = 0 0 0.0
> pull_rate1 = 0
> pull_k1 = <<<force-constant>>>
> pull_vec1 = 0 0 0
>
> -- original message --
>
>
> The manual does discuss restraining to a plane, but this must be the plane
> in which the atom is already present.
>
> [ position_restraints ]
> ; ai funct fc
> ...
> 3 1 1000 0 0
>
> How about restraining the atom to some other plane? For example, how about
> restraining a phosphate group initially at the lipid-water interface to the
> bilayer center (for whatever fancy reasons) ? Won't this require pulling
> the
> atom to that plane first?
>
> If it can be achieved using the position restraints alone, it is not clear
> to me how to do this?
>
> -Maria
>
>
>
> On Tue, Sep 7, 2010 at 12:53 PM, Justin A. Lemkul <jalemkul at vt.edu>
> wrote:
>
> maria goranovic wrote:
>
> I want to restrain certain atoms of my simulation to the plane
> perpendicular to the bilayer normal, and at the bilayer center. Can
> someone
> please provide a quick guide on how to do this? I read the pull-code
> options, but restraining to a plane did not seem possible?
>
>
> You can restrain atoms to a plane using position restraints. See the
> manual for an example.
> -Justin
> --
> Maria G.
> Technical University of Denmark
> Copenhagen
>
> --
> gmx-users mailing list gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/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/Support/Mailing_Lists
>
--
Maria G.
Technical University of Denmark
Copenhagen
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20100909/81046af4/attachment.html>
More information about the gromacs.org_gmx-users
mailing list