[gmx-users] Simulated Annealing
Martina Bertsch, PhD
mbe404 at lulu.it.northwestern.edu
Mon Jan 19 08:30:02 CET 2004
I would like to perform the following simulated annealing protocol on a
GPCR model in vacuum, while keeping the C-alpha atoms of the
transmembrane (TM) helices restrained (fixed?):
1) Quick heating from 300-600 K in 1 ps with a 1 fs step.
2) Equilibration at 600 K for 100 ps.
3) Slow cooling to 500 K.
4) Equilibration at 500 K for 200 ps.
5) Slow cooling to 400 K.
6) Equilibration at 400 K for 200 ps.
7) Slow cooling to 300 K.
8) Equilibration at 300 K for 200 ps with gradual release of the
restrained (fixed?) TM C-alpha atoms.
9) Minimization of the resulting structure.
My questions are:
a) Is there a way to program the whole protocol in a single mdp? Which
mdp options do I use?
In the past, I used a simple SA cooling from 600 to 0 K in 1 ns (see
enclosed mdp), but the papers suggest that a more gradual cooling is
2) During this process, I would like to keep just the transmembrane
C-alpha atoms fixed, e.g., C-alpha in residues 6-30, 43-62, etc. Should I:
make_ndx -f GPCR_postsd.gro -o index.ndx
> r 6-30 & a CA
> r 43-62 & a CA
I seem to recall reading that an atom can only belong to one group. Does
that mean that I have to split the original [ C-alpha ] into the [
and the [nonTM_C-alpha] and delete the original [ C-alpha ] group?
Wouldn't that also imply that I should get rid of the [ backbone ] group
-- that sounds a bit suspicious?
3) Once I obtain the [ TM_C-alpha ] group, how do I fix it in the mdp?
Should I just define:
freeze_dim y y y
Or, should I prepare a restraint file TM_posres.itp and define
-DPOSRES_TM in the mdp?
4) How would I encode the gradual release of the [TM_C-alpha] group at
the last stage of the protocol?
5) If I did not want to fix any parts of the structure, how would I
obtain the energy and rmsd for parts of the structure (with respect to
the beginning structure), e.g., for TM helix 5, TM helix 7, etc?
I know that I could define groups TM1-TM7 in the index.ndx, and then for
the rmsd calculation:
g_rms -s gpcr_SA.tpr -f gpcr_SA.trr -o TM7.xvg
and select the desired group at the prompt.
I am not sure how to obtain energy profiles for certain groups.
The "simple" SA protocol:
title = MT1 Simulated Annealing
cpp = /lib/cpp
constraints = all-bonds
integrator = md
dt = 0.001 ; ps !
nsteps = 1000000 ; total 1 ns.
;nstcomm = 1
nstxout = 2000
nstvout = 2000
nstfout = 2000
nstlog = 2000
nstenergy = 2000
nstxtcout = 2000
nstlist = 10
ns_type = grid
rlist = 1.0
coulombtype = PME
rcoulomb = 1.0
rvdw = 1.0
pbc = xyz
comm_mode = linear
nstcomm = 1
fourierspacing = 0.15
pme_order = 4
optimize_fft = yes
comm_grps = protein
; Berendsen temperature coupling is on in one group
Tcoupl = berendsen
tau_t = 0.1
tc-grps = protein
ref_t = 600
; Energy monitoring
energygrps = protein
; Pressure coupling is not on
Pcoupl = no
pcoupltype = isotropic
tau_p = 5
compressibility = 4.5e-5
ref_p = 1.0
; Simulated annealing
annealing = yes
zero_temp_time = 1000000 ; ps
; Generate velocites is on at 600 K.
gen_vel = yes
gen_temp = 600.0
gen_seed = 280572
More information about the gromacs.org_gmx-users