[gmx-users] Basic GROMACS questions

[Justin Lemkul]

On 12/8/14 10:47 AM, Simone Bolognini wrote:
> Thank you Justin. Your answer was very helpful. Since I tries like 3-4
> GROMACS tutorials I wanted to try to simulate dynamics with my own
> parameters. I'm trying to simulate the chicken prion protein fragment 1U3M
> using gromos53a6 as force field. I followed the typical protocol for
> dynamics (vacuum minimization, addition of solvent and neutralization, in
> solvent minimization, position restrained dynamics and production). Here
> are my .mdp files for the various tasks:
>
> ;VACUUM MINIMIZATION
> define = -DFLEXIBLE
> integrator         = cg
> emtol = 50
> nstcgsteep = 1000
> nsteps = -1
>
> nstenergy = 1
> energygrps = System
>
> nstlist = 1
> nstype = simple
> pbc = no
>
> coulombtype = cut-off
> vdwtype = cut-off
> rvdw = 1
> rlist = 1
> rcoulomb = 1
>
> actually I could have used PME electrostatics with pbc here
> emtol is set in such a way that convergence is really reached (I tried with
> different values and looked at the convergence of the potential energy)
>

I don't know why you want to have PME and PBC in vacuum.  There's no point.  In 
vacuo minimization can also be dangerous in some cases, when the identified 
minimum has become unphysical...

>
> ;SOLVENT MINIMIZATION
> define = -DFLEXIBLE
>
> integrator = cg
> emtol = 500
> nstcgsteep = 1000
> nsteps = -1
>
> nstenergy = 1
> energygrps = System
>
> nstlist = 10
> nstype = grid
> pbc = xyz
>
> coulombtype = PME
> vdwtype = cut-off
> rvdw = 1
> rlist = 1
> rcoulomb = 1
> fourierspacing = 0.15
> pme-order = 4
> ewald-rtol = 1e-5
> optimize-fft = yes
>
>
> ;PRD
> define = -DPOSRES
> constraints     = all-bonds
> constraint_algorithm = LINCS
>
> integrator = md-vv-avek
> dt = 0.002
> nsteps = 50000
> nstcomm = 10
>
> nstenergy = 10
> nstxout = 500
> nstxtcout      = 500
> nstvout = 5000
> nstfout = 0
> nstlog = 10
>
> nstlist = 10
> nstype = grid
> pbc = xyz
>
> coulombtype = PME
> vdwtype = cut-off
> rvdw = 1.
> rlist = 1.
> rcoulomb = 1.
> fourierspacing = 0.16
> pme-order = 4
> ewald-rtol = 1e-5
> optimize-fft = yes
>
> tcoupl       = Berendsen
> tau_t       = 0.5  0.5
> tc_grps       = protein non-protein
> ref_t       = 300 300
> nsttcouple = 1
>
> pcoupl = Berendsen
> pcoupltype = isotropic
> tau_p = 1.0
> compressibility = 4.5e-5
> ref_p = 1.0
> refcoord_scaling = all
> nstpcouple = 1
>
> gen_vel = yes
> gen_temp = 300
> gen_seed = -1
>
>
> ;PRODUCTION
> constraints     = all-bonds
> constraint_algorithm = lincs
> lincs_iter = 1
> lincs_order = 4
>
> integrator = md-vv-avek
> dt       = 0.002
> nsteps = 500000
> nstcomm = 10
>
> nstenergy = 10
> nstxout = 1000
> nstxtcout = 1000
> nstvout = 1000
> nstlog = 1000
>
> nstlist = 5
> nstype = grid
> pbc = xyz
>
> coulombtype = PME
> vdwtype       = cut-off
> rvdw = 1.0
> rlist = 1.0
> rcoulomb = 1.0
> fourierspacing = 0.16
> pme-order = 4
> ewald-rtol = 1e-5
> optimize-fft = yes
>
> tcoupl       = nose-hoover  ; nose-hoover extended ensamble
> tau_t       = 0.5 0.5   ; relaxation-time
> tc_grps       = protein non-protein
> ref_t       = 300 300
>
> pcoupl       = Parrinello-Rahman
> pcoupltype      = isotropic
> tau_p       = 2.0
> compressibility = 4.5e-5
> ref_p       = 1.0
> refcoord_scaling = all
> gen_vel      = no
>
>
> Thank you again for your help!!!
>

Note that a plain 1.0 nm cutoff is incorrect for Gromos96 force fields.  A 
0.9/1.4 twin-range setup is how the force field was parametrized.

-Justin

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

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

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