[gmx-users] How to define pull-init for two pull groups
Rini Gupta
rinisgupta at gmail.com
Wed Jan 15 02:23:42 CET 2014
Dear gmx-users,
Thanks for the nice reply.
I have corrected the pull code as:
; Pull code
pull = umbrella
pull_geometry = position ;
pull_dim = N N Y
pull_start = yes ; define initial COM distance > 0
pull_ngroups = 2
pull_group0 = cbilayer ; center of carbon's of bilayer
pull_group1 = phosup ; phosphate atom of 16th PEPC in top
pull_pbcatom1 = 0
pull_vec1 = 0.0 0.0 -1.0
pull_init1 = 0 0 0
pull_rate1 = 0.01 ; 0.01 nm per ps = 10 nm per ns
pull_k1 = 1000 ; kJ mol^-1 nm^-2
pull_group2 = phoslow ; phosphate atom of 56th PEPC in bottom
pull_pbcatom2 = 0
pull_vec2 = 0.0 0.0 1.0
pull_init2 = 0 0 0
pull_rate2 = 0.01 ; 0.01 nm per ps = 10 nm per ns
pull_k2 = 1000 ; kJ mol^-1 nm^-2
Now, after extracting trajectory frames and measuring the distances between
COM distances between two pull groups and a reference group, I found that
values of summary_distances.dat are first decreasing through roughly half
of the frames and then start increasing till final frame.
The snapshots are also matching the output.
If this is right for a lipid crossing through bilayer then in order to
generate starting configurations at windows spacing of 0.06 nm I ran the
script (setup-umbrella.py) link provided in tutorial.
But, the output file (caught_output.txt) which contain the list of initial
frames with distances
select only frames after half of the times frames.
i.e If I provide a .dat file corresponding to 400 time frames it produces
the list of frames at a desired spacing excluding the frames where distance
between reference and pull group is decreasing . I got :
Creating frame-specific output for files:
run-umbrella.sh
frame dist d_dist
0 3.082 NA
6 3.137 0.055
222 3.215 0.078
224 3.272 0.057
228 3.320 0.049
230 3.380 0.060
231 3.427 0.047
237 3.495 0.068
239 3.557 0.062
240 3.619 0.062
242 3.720 0.102
257 3.781 0.061
271 3.835 0.053
272 3.894 0.059
283 3.935 0.041
286 4.010 0.075
288 4.063 0.053
299 4.113 0.049
300 4.188 0.075
302 4.243 0.055
305 4.303 0.060
311 4.397 0.094
312 4.431 0.034
If this is case, how it will generate a sufficient overlap between adjacent
windows ?
Also, script works only for single lipid pulling out of a layer. It is
possible to use it for two lipids simultaneously crossing the bilayer in
different directions.
Thanks and Regards,
Rini
<gmx-users at gromacs.org>
>
>
>
>
>
> On 1/13/14, 7:40 PM, Rini Gupta wrote:
>
>> Dear gmx-users,
>>
>> I want to calculate Potential of Mean Force (PMF) of a lipid in a bilayer
>> by umbrella sampling method using coarse grain molecular dynamics
>> simulations. My simulation system consists of
>> 76 (PEPC and GDPE lipids mixture), 2495 water molecules and 16 other
>> molecules which are acting as a surface attached to this bilayer system.
>>
>>
>> In order to generate starting structures for the pulling simulations,
>> after
>> minimization, I have equilibrated the system for 25 ns under NVT ensemble
>> followed by 25 ns NPT run with position restraints applied to both lipids
>> and surface molecules. Following this, I have ran production simulations
>> under NPT ensemble for another 25 ns. During production phase, position
>> restraints were removed from all lipids except surface molecules.
>> I want to calculate PMF for pulling out two lipids simultaneously one
>> in
>> top and one in bottom layer of this system.
>> Following Justin tutorial, I understand that I need to define two pull
>> groups and one reference group for my problem. For this, I have chosen all
>> carbon's of lipid bilayer as reference group and Phosphate atom of two
>> lipids which are placed in top and bottom layers respectively,
>> as pull group1=phosup and pull group2 =phoslow
>> My problem is how to define pull-init for two pull groups:
>> I have to directly place the initial coordinates of two selected lipids as
>> pull-init or first I need to assume that bilayer center is shifted to
>> (0,0,0 ) and then put initial cordinates as pull-init.
>> Please let me know if anything is wrong here.
>> Following second approach, I have used the following grompp.mdp file for
>> pulling simulation:
>>
>> ; Pull code
>> pull = umbrella
>> pull_geometry = position ;
>> pull_dim = N N Y
>> pull_start = yes ; define initial COM distance > 0
>> pull_ngroups = 2
>> pull_group0 = cbilayer ; center of carbon's of bilayer
>> pull_group1 = phosup ; phosphate atom of 16th PEPC in top
>> pull_pbcatom1 = 0
>> pull_vec1 = 0.0 0.0 -1.0
>> pull_init1 = 0.058 -1.88226 -1.91169
>> pull_rate1 = 0.01 ; 0.01 nm per ps = 10 nm per ns
>> pull_k1 = 1000 ; kJ mol^-1 nm^-2
>> pull_group2 = phoslow ; phosphate atom of 56th PEPC in bottom
>> pull_pbcatom2 = 0
>> pull_vec2 = 0.0 0.0 1.0
>> pull_init2 = 1.778 1.20774 2.20831
>> pull_rate2 = 0.01 ; 0.01 nm per ps = 10 nm per ns
>> pull_k2 = 1000 ; kJ mol^-1 nm^-2
>>
>>
> One problem is that you are using both "pull_start = yes" and specifying
> non-zero values for pull_init[12]. What you are doing is specifying the
> reference distance as the initial COM distance + pull_init[12]. If you are
> using "pull_start = yes," set pull_init[12] to a zero vector. If you
> manually specify the actual COM distance in pull_init[12], then set
> "pull_start = no." The trajectories indicate that you are probably getting
> something close to what you want to achieve, but likely the forces and
> resulting output distances are going to be somewhat screwy.
>
>
>
>> After this, following tutorial , I extracted the separate trajectory
>> frames
>> using trjconv and then in order to measure the COM distances between two
>> pull groups and a reference group
>> I have created two groups.txt files and as a result I got two
>> summary_distances.dat files corresponding to two pulling lipids with
>> respect to a reference group.
>> However, looking at the contents of both files I do not observe the any
>> definite decreasing or increasing COM distances between pull group and
>> reference.
>> But, on observing the snapshots,
>> I found that pull groups are moving in right direction.
>> If I am pulling a lipid across bilayer,
>> these COM distance between pull group and reference group should first
>> decrease and then increase right?
>> Please tell me what I am doing wrong here.
>>
>> How to take some of these configurations as starting window structures
>> for
>> next step of pulling simulations?
>>
>>
> That depends on what you define as a reaction coordinate. If you are
> interested in sampling states during which both lipids are traversing the
> bilayer, then you need configurations that are appropriate for both lipids
> along the reaction coordinate. There's nothing really special about it,
> except that you likely have two target distances to simultaneously satisfy.
>
> -Justin
>
> --
> ==================================================
>
> Justin A. Lemkul, Ph.D.
> Postdoctoral Fellow
>
> Department of Pharmaceutical Sciences
> School of Pharmacy
> Health Sciences Facility II, Room 601
> University of Maryland, Baltimore
> 20 Penn St.
> Baltimore, MD 21201
>
> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>
> ==================================================
> --
> Gromacs Users mailing list
>
> * Please search the archive at http://www.gromacs.org/
> Support/Mailing_Lists/GMX-Users_List before posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>
>
More information about the gromacs.org_gmx-users
mailing list