[gmx-users] Umbrella Sampling on Lipid Bilayers, Pulling Issue and Umbrella Sampling Window Simulation Issues.

Justin Lemkul jalemkul at vt.edu
Mon Nov 4 20:43:04 CET 2019



On 11/1/19 2:09 PM, Alexander Phillips wrote:
> Hello Everyone,
>
> I am currently running a CG MARTINI simulation on a bilayer of POPC Lipid
> in Water. I also have a single vitamin E in the mixture aswell (inside of
> the POPC bilayer). I am trying to perform umbrella sampling by pulling the
> vitamin E and measuring it's binding energy to POPC. I am basing most of my
> work off of the tutorial made by Justin Lemkul. The plan for my umbrella
> sampling is to pull the vitamin E vertically out of the membrane  (from the
> center of the membrane to some distance into the water). However I am
> having issues with performing umbrella sampling. I believe the issue is one
> of two things (or perhaps a mixture of both):
> 1) The simulation where I pull the Vitamin E out of POPC is happenning too
> quickly and I cannot gather enough frames to get valid input windows for
> the umbrella sampling portion.
> 2) The umbrella sampling portion is not keeping the Vitamin E around the
> starting portion and it falls back into the membrane.
>
> I will attach the mdp files I used for both, and hopefully someone can give
> me some advice about how to do this / fix it.

It may be a combination of things. Your pull rate is quite fast and your 
force constant is very small. I suggest experimenting with slower 
pulling and/or a stiffer spring constant.

> I am running in Gromacs 4.6.5 so a few things are changed from the files
> given in Justin's tutorial, however I believe ultimately it's about the
> same.

Any reason to use such antiquated software? The modern versions will be 
MUCH faster (and are actually supported by the developers in case of 
problems).

-Justin

> Here is the mdp for the Pulling:
>
>
> integrator         = md
> dt                        = 0.03
> nsteps               = 250000
> nstcomm          = 10
>
> comm-grps      =
> nstxout             = 5000
> nstvout             = 5000
> nstfout              = 5000
> nstlog                = 1000  ; Output frequency for energies to log file
> nstenergy         = 100   ; Output frequency for energies to energy file
> nstxtcout          = 1000  ; Output frequency for .xtc file
> xtc_precision   = 5
> xtc-grps             =
>
> energygrps       = Lipids Solvent PULLING_ATOC
> nstlist                  = 10
> ns_type               = grid
> pbc                      = xyz
> rlist                       = 1.4
> coulombtype    = Shift  ;Reaction_field (for use with Verlet-pairlist) ;PME
> (especially with polarizable water)
> rcoulomb_switch          = 0.0
> rcoulomb           = 1.2
> epsilon_r            = 15   ; 2.5 (with polarizable water)
> vdw_type            = Shift  ;cutoff (for use with Verlet-pairlist)
> rvdw_switch      = 0.9
> rvdw                    = 1.2  ;1.1 (for use with Verlet-pairlist)
>
> tcoupl                 = v-rescale
> tc-grps                = Lipids Solvent PULLING_ATOC
> tau_t                    = 1.0  1.0 1.0
> ref_t                     = 310 310 310
> Pcoupl                = parrinello-rahman
> Pcoupltype        = semiisotropic
> tau_p                   = 12.0 12.0 12.0 ;parrinello-rahman is more stable
> with larger tau-p, DdJ, 20130422
> compressibility = 3e-4 3e-4 3e-4
> ref_p                    = 1.0 1.0 1.0
>
> gen_vel               = no
> gen_temp           = 310
> gen_seed            = 473529
>
> constraints         = none
> constraint_algorithm  = Lincs
> continuation      = no
> lincs_order          = 4
> lincs_warnangle  = 30
>
> ; Pull Code
> pull = umbrella
> pull-geometry = direction-periodic
> pull-start = yes
> pull-nstxout = 10
> pull-nstfout = 10
> pull-ngroups = 1
> pull-group0 = Lipids
> pull-group1 = PULLING_ATOC
> pull-vec1 = 0 0 1
> pull-rate1 = 0.01
> pull-k1 = 20
>
>
>
>
> ---------------------------------------------------------------- End of MDP
> 1 -----------------------
>
> The second MDP file (the one used for the Umbrella Sampling production run)
> is listed here
>
> integrator             = md
> dt                            = 0.002
> nsteps                   = 10000000
> nstcomm              = 10
> comm-grps           =
>
> nstxout                  = 50000
> nstvout                  = 50000
> nstfout                   = 5000
> nstlog                     = 1000  ; Output frequency for energies to log
> file
> nstenergy              = 5000  ; Output frequency for energies to energy
> file
> nstxtcout               = 5000  ; Output frequency for .xtc file
> xtc_precision        = 100
> xtc-grps                 =
> energygrps           = Lipids Solvent PULLING_ATOC
>
> nstlist                = 5
> ns_type             = grid
> pbc                    = xyz
> rlist                    = 1.4
>
> coulombtype       = Shift  ;Reaction_field (for use with Verlet-pairlist)
> ;PME (especially with polarizable water)
> rcoulomb_switch = 0.0
> rcoulomb              = 1.2
> epsilon_r               = 15   ; 2.5 (with polarizable water)
> vdw_type              = Shift  ;cutoff (for use with Verlet-pairlist)
> rvdw_switch         = 0.9
> rvdw                       = 1.2  ;1.1 (for use with Verlet-pairlist)
>
> tcoupl                    = v-rescale
> tc-grps                   = Lipids Solvent PULLING_ATOC
> tau_t                       = 1.0  1.0 1.0
> ref_t                        = 310 310 310
> Pcoupl                   = parrinello-rahman
> Pcoupltype           = semiisotropic
> tau_p                      = 12.0 12.0 12.0 ;parrinello-rahman is more
> stable with larger tau-p, DdJ, 20130422
> compressibility    = 3e-4 3e-4 3e-4
> ref_p                       = 1.0 1.0 1.0
>
> gen_vel                      = no
> gen_temp                 = 310
> gen_seed                  = 473529
>
> constraints                   = none
> constraint_algorithm = Lincs
> continuation                = no
> lincs_order                    = 4
> lincs_warnangle          = 30
>
> ; Pull Code
> pull                     = umbrella
> pull-geometry = direction-periodic
> pull-start           = yes
> pull-nstxout     = 10
> pull-nstfout      = 10
> pull-ngroups    = 1
> pull-group0      = Lipids
> pull-group1      = PULLING_ATOC
> pull-vec1            = 0 0 1
> pull-rate1           = 0
> pull-k1                 = 1
>
> ---------------------------- END of MDP 2 -------------------
>
>
> If you have any questions, please ask and I will try my best to answer
> them, Thank you,
> Alex

-- 
==================================================

Justin A. Lemkul, Ph.D.
Assistant Professor
Office: 301 Fralin Hall
Lab: 303 Engel Hall

Virginia Tech Department of Biochemistry
340 West Campus Dr.
Blacksburg, VA 24061

jalemkul at vt.edu | (540) 231-3129
http://www.thelemkullab.com

==================================================



More information about the gromacs.org_gmx-users mailing list