[gmx-users] To get PMF using pull geometry of cylinder
Lin, Dejun
Dejun_Lin at URMC.Rochester.edu
Thu May 26 18:21:43 CEST 2011
Ah...Sorry if I confused you. Actually what I meant is if I set the pull code, e.g., as:
pull_geometry = position
pull_dim = N N Y
pull_rate1 = 0.1
pull_init1 = 3.0
pull_vec1 = 0 0 1
(where pull rate is NOT zero)
Then the position = pull_init + time*pull_rate*pull_vec1 and since pull_vec1 is z-component-positive (0, 0, 1), the position between the reference and the pull group increases monotonically and that's not what I want. I want the distance to decrease so the pull_vec1 should be (0,0,-1). Am I right?
Thanks,
Dejun
________________________________________
From: gmx-users-bounces at gromacs.org [gmx-users-bounces at gromacs.org] On Behalf Of Justin A. Lemkul [jalemkul at vt.edu]
Sent: Thursday, May 26, 2011 11:26 AM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] To get PMF using pull geometry of cylinder
Lin, Dejun wrote:
> Thank you Justin. But I'm still confused about the pull_vec1. Is the vector
> pointing from the reference or from the pull group for all the geometry
> settings? For instance, if I use pull_geometry = position, pull_dim = N N Y,
> pull_rate1 = 0.1, pull_init1 = 3.0, pull_vec1 = 0 0 1, the position =
> pull_init + time*pull_rate*pull_vec1 according to the manual and it's
> actually pulling the pull group away from the reference because the vector
> has a positive z component.
>
The reference group serves as the reference :) Thus, the vector connecting your
reference and pull group is (0,0,1) since the pull group has a larger
z-coordinate. You're right about the position equation, but you've got a
pull_rate1 of zero, so the added term drops out and position = pull_init1.
> BTW, I did consider using position geometry but can I get the PMF using
> pull_rate NOT equal to zero? The reason I turned to cylinder is that with
> "distance" geometry the membrane concaved as the micelle was restrained to a
> very short distance (any distance below the critical distance where the
> fusion was supposed to occur) without fusing and this would persist for some
> hundred nanoseconds.
>
To obtain the PMF, the pull rate should be zero, using several (many)
independent simulations along the desired reaction coordinate. See the tutorial:
http://www.gromacs.org/Documentation/Tutorials#Pull_Code_and_Umbrella_Sampling
-Justin
> Thanks again for your help! Dejun
>
> ________________________________________ From: gmx-users-bounces at gromacs.org
> [gmx-users-bounces at gromacs.org] On Behalf Of Justin A. Lemkul
> [jalemkul at vt.edu] Sent: Wednesday, May 25, 2011 7:44 PM To: Discussion list
> for GROMACS users Subject: Re: [gmx-users] To get PMF using pull geometry of
> cylinder
>
> Lin, Dejun wrote:
>> Hi all,
>>
>> I want to do an umbrella sampling to calculate the potential of mean force
>> of a micelle fusing to a lipid membrane. In the starting configuration, the
>> membrane normal is set in z direction and the membrane is under the
>> micelle (z(COM-membrane) < z(COM-micelle)). At the beginning, I used
>> pull_geometry = distance but the sampling took to long to converge. And
>> then I considered using cylinder instead and set up the pull code as: pull
>> = umbrella pull_geometry = cylinder pull_dim = N N Y pull_start =
>> no pull_ngroups = 1 pull_group0 = membrane pull_group1 = micelle
>> pull_vec1 = 0.0 0.0 -1.0 pull_r1 = 2.3 pull_r0 = 2.5 pull_init1
>> = 1.8 pull_rate1 = 0.0 pull_k1 = 1000 pull_nstxout = 1000
>> pull_nstfout = 1000
>>
>> And the output files (px.xvg and pf.xvg) look like: (px.xvg) 0.0000
>> 23.3397 4.90267 100.0000 0.742426 1.86002
>> 200.0000 0.737842 1.86232 300.0000
>> 0.762641 1.83985 400.0000 21.6814
>> 1.86777 500.0000 0.771017 1.82482 600.0000
>> 0.789036 1.83254 700.0000 0.793745
>> 1.83503 800.0000 0.804732 1.84217 900.0000
>> 0.80137 1.84166 1000.0000 0.817979 1.83546
>> 1100.0000 0.806272 1.83535 1200.0000 0.792355
>> 1.80372 1300.0000 0.786839 1.8438 1400.0000 0.798686
>> 1.82625 1500.0000 0.802859 1.84138 1600.0000 0.803844
>> 1.84379
>>
>> (pf.xvg) 100.0000 -480.158 200.0000 -498.547
>> 300.0000 -318.839 400.0000 -542.14 500.0000
>> -198.526 600.0000 -260.291 700.0000 -280.22
>> 800.0000 -337.385 900.0000 -333.317 1000.0000 -283.644
>> 1100.0000 -282.761 1200.0000 -29.7684 1300.0000 -350.409 1400.0000
>> -209.991 1500.0000 -331.029 1600.0000 -350.309
>>
>> And I visualized the trajectory in VMD and found that an "explosion" of the
>> micelle (with all its residues scattered in the box). I then changed the
>> pull_vec1 = 0.0 0.0 -1.0 to pull_vec1 = 0.0 0.0 1.0 and the "explosion"
>> was gone and everything seemed to work. But I'm not quite sure if I'm
>> setting up the pull-code right because if z(COM-membrane) < z(COM-micelle),
>> then pull_vec1 = 0.0 0.0 -1.0 should be the right choice. Basically, I
>> think my
>
> The negative is incorrect. If the micelle is "above" the membrane (i.e., it
> has a greater z-coordinate), then the reference vector should be positive
> with respect to the reference position. It would appear from the pullx.xvg
> file that you've got some weird PBC crossing going on, since the reference
> positions oscillate dramatically at the beginning of the simulation.
>
>> pull-code is restraining the COM of the micelle to the COM of a "cylinder"
>> of the membrane at a distance of 1.8nm, with the vector pointing from
>> COM(micelle) to COM(cylinder) being (0,0,-1). So can anyone explain what
>> caused the explosion or am I setting up the pull-code right? And what are
>> the 2nd and 3rd column in the px.xvg output in this case?
>>
>
> For such a simple restraint, I don't think there is any need for cylinder
> geometry. You're not using different weighting. Any of the other
> pull_geometry settings should work, with "distance" probably being the most
> intuitive. I suspect it would be very easy to use the "position" setting as
> well.
>
> The output in the pullx.xvg file should be labeled in the header. The
> generic output is (x,y,z) of the reference group, followed by delta(x,y,z)
> for the restrained group, i.e. how far away the restrained group is along
> each coordinate axis. If you're only restraining/pulling along the z-axis,
> then the only terms written are z and dz, since, in principle, x, y, dx, and
> dy are all fixed, with the last two terms being zero.
>
> -Justin
>
> -- ========================================
>
> Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee
> Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu |
> (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
>
> ======================================== -- 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 the www interface or send
> it to gmx-users-request at gromacs.org. Can't post? Read
> http://www.gromacs.org/Support/Mailing_Lists
>
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
--
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 the
www interface or send it to gmx-users-request at gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
More information about the gromacs.org_gmx-users
mailing list