On 6/27/12 9:36 AM, Steven Neumann wrote:
On Wed, Jun 27, 2012 at 1:51 PM, Justin A. Lemkul <jalem...@vt.edu> wrote:
On 6/27/12 7:48 AM, Steven Neumann wrote:
Dear Gmx Users,
I obtained a protein-ligand complex from 100ns simulation. Now I am
pulling my ligand away from the protein after the energy minimzation
in water and equilibration of 100ps (two coupling baths: Protein,
LIG_Water_and_ions).
Then I proceed my pulling :
grompp -f pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o
pull.tpr
mdrun -s pull.tpr -deffnm pull
title = Umbrella pulling simulation
define = -DPOSRES
; Run parameters
integrator = md
dt = 0.002
tinit = 0
nsteps = 500000 ; 1 ns
nstcomm = 10
; Output parameters
nstxout = 0
nstvout = 0
nstfout = 500
nstxtcout = 1000 ; every 1 ps
nstenergy = 500
; Bond parameters
constraint_algorithm = lincs
constraints = all-bonds
continuation = yes ; continuing from NPT
; Single-range cutoff scheme
nstlist = 5
ns_type = grid
rlist = 1.4
rcoulomb = 1.4
rvdw = 1.2
vdwtype = Switch
rvdw-switch = 1.0
; PME electrostatics parameters
coulombtype = PME
fourierspacing = 0.12
fourier_nx = 0
fourier_ny = 0
fourier_nz = 0
pme_order = 4
ewald_rtol = 1e-5
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 = 298 298 ; reference
temperature, one for each group, in K
; Pressure coupling is on
Pcoupl = Parrinello-Rahman
pcoupltype = isotropic
tau_p = 2.0
compressibility = 4.5e-5
ref_p = 1.0
; Generate velocities is off
gen_vel = no
; Periodic boundary conditions are on in all directions
pbc = xyz
; Long-range dispersion correction
DispCorr = EnerPres
; Pull code
pull = umbrella
pull_geometry = distance ; simple distance increase
pull_dim = N N Y
pull_start = yes ; define initial COM distance > 0
pull_ngroups = 1
pull_group0 = Protein
pull_group1 = LIG
pull_rate1 = 0.004 ; 0.004 nm per ps = 4 nm per ns
pull_k1 = 500 ; kJ mol^-1 nm^-2
I run 3 pulling simulations with the same mdp and I obtain 3
different profiles (Force vs time). Then I used 2xlonger pulling with
the same pulling distance and I run 3 simulations again. Each time I
obtain different profile. Can anyone explain me this? I am using
velocities from npt simulation as above (gen_vel = no and continuation
= yes) so I presume the output should be similar. Please, advice.
I assume you're passing a checkpoint file to grompp? If you're relying on
velocities from the .gro file, they are of insufficient precision to
guarantee proper continuation.
Thank you Justin. I am using according to your tutorial:
grompp -f pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o pull.tpr
mdrun -s pull.tpr -deffnm pull
Would you suggest:
grompp -f pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o pull.tpr
mdrun -s pull.tpr -cpi npt.cpt -deffnm pull ??
No, I would not, especially if the NPT run uses position restraints, in which
case the two phases are different. I missed the command line in the earlier
message. What you are doing makes sense.
Profiles do not vary slightly - the maximum pulling force (breaking
point) varies from 290 to 500 kJ/mol nm2 which is really a lot.
Consult the points below and watch your trajectories. If you're getting
different forces, your ligands are experiencing different interactions. SMD is
a path-dependent, non-equilibrium process. Good sampling and a justifiable path
are key.
-Justin
Small variations are inherent in any simulation set, and in the case of
pulling, small changes (though intentional) are the basis for Jarzynski's
method. In any case, all MD simulations are chaotic and so it depends on
what your definition of "different" is in the context of whether or not
there are meaningful changes imparted through the course of each simulation.
Also note that in the absence of the -reprod flag, the same .tpr file may
result in a slightly different outcome. The implications of these outcomes
are limited by sampling; the ensemble should converge with sufficient time
and/or replicates. For non-equilibrium processes like pulling, convergence
is probably harder, but again you have to ask whether the differences are
meaningful.
http://www.gromacs.org/Documentation/Terminology/Reproducibility
-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
--
========================================
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
