On 8/21/12 11:09 AM, Steven Neumann wrote:
On Tue, Aug 21, 2012 at 3:48 PM, Justin Lemkul <jalem...@vt.edu> wrote:
On 8/21/12 10:42 AM, Steven Neumann wrote:
Please see the example of the plot from US simulations and WHAM:
http://speedy.sh/Ecr3A/PMF.JPG
First grompp of frame 0 corresponds to 0.8 nm - this is what was shown
by grompp at the end.
The mdp file:
; Run parameters
define = -DPOSRES_T
integrator = md ; leap-frog integrator
nsteps = 25000000 ; 100ns
dt = 0.002 ; 2 fs
tinit = 0
nstcomm = 10
; Output control
nstxout = 0 ; save coordinates every 100 ps
nstvout = 0 ; save velocities every
nstxtcout = 50000 ; every 10 ps
nstenergy = 1000
; Bond parameters
continuation = no ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = all-bonds ; all bonds (even heavy atom-H bonds)
constrained
; Neighborsearching
ns_type = grid ; search neighboring grid cells
nstlist = 5 ; 10 fs
vdwtype = Switch
rvdw-switch = 1.0
rlist = 1.4 ; short-range neighborlist cutoff (in nm)
rcoulomb = 1.4 ; short-range electrostatic cutoff (in nm)
rvdw = 1.2 ; short-range van der Waals cutoff (in nm)
ewald_rtol = 1e-5 ; relative strenght of the Ewald-shifted
potential rcoulomb
; Electrostatics
coulombtype = PME ; Particle Mesh Ewald for long-range
electrostatics
pme_order = 4 ; cubic interpolation
fourierspacing = 0.12 ; grid spacing for FFT
fourier_nx = 0
fourier_ny = 0
fourier_nz = 0
optimize_fft = yes
; Temperature coupling is on
tcoupl = V-rescale ; modified Berendsen
thermostat
tc_grps = Protein LIG_Water_and_ions ; two coupling groups - more
accurate
tau_t = 0.1 0.1 ; time constant, in ps
ref_t = 318 318 ; reference temperature,
one for each group, in K
; Pressure coupling is on
pcoupl = Parrinello-Rahman ; pressure coupling is on for
NPT
pcoupltype = isotropic ; uniform scaling of box
vectors
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 ; reference pressure, in bar
compressibility = 4.5e-5 ; isothermal
compressibility of water, bar^-1
; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes ; assign velocities from Maxwell distribution
gen_temp = 318 ; temperature for Maxwell distribution
gen_seed = -1 ; generate a random seed
; These options remove COM motion of the system
; Pull code
pull = umbrella
pull_geometry = distance
pull_dim = N N Y
pull_start = yes
pull_ngroups = 1
pull_group0 = Protein
pull_group1 = LIG
pull_init1 = 0
pull_rate1 = 0.0
pull_k1 = 500 ; kJ mol^-1 nm^-2
pull_nstxout = 1000 ; every 2 ps
pull_nstfout = 1000 ; every 2 ps
Based on these settings you're allowing grompp to set the reference distance
to whatever it finds in the coordinate file. It seems clear to me that the
sampling indicates what I said before - you have an energy minimum somewhere
other than where you "started" with. What that state corresponds to
relative to what you think is going on is for you to decide based on the
nature of your system. Whatever is occurring at 0.6 nm of COM separation is
of particular interest, since the energy minimum is so distinct.
So based on this the deltaG will correspond to -5.22 as the initial
state was taken at 0.4 nm corresponding to 0 kcal/mol as the moment
corresponding to the minimum is the coordinate from SMD where last
hydrogen bond was broken. Would you agree?
Based on the very little information I have, no. It would appear that the 0.4
nm separation is in fact some metastable state and the true energy minimum is at
0.6 nm of COM separation. What's going on at that location?
I hope you're doing a thorough analysis of convergence if you're generating
velocities at the outset of each run, and removing unequilibrated frames
from your analysis.
When I use WHAM I skip first 200 ps of each window as the equilibration.
That seems fairly short, especially given the generation of velocities in
conjunction with the Parrinello-Rahman barostat, which can be very temperamental.
-Justin
--
========================================
Justin A. Lemkul, Ph.D.
Research Scientist
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@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Only plain text messages are allowed!
* 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-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
