[gmx-users] topology and parameter set up
Justin Lemkul
jalemkul at vt.edu
Mon Feb 2 14:40:41 CET 2015
On 2/2/15 8:20 AM, Jennifer Vo wrote:
> Dear Justin,
> Many thanks for your help. My em.mdp is
> integrator = steep
> emtol = 100.0
> emstep = 0.01
> nsteps = 100000
> coulombtype = PME
> pbc = xyz
>
> My nvt.mdp is
> title = Protein-ligand complex NVT equilibration
> define = -DPOSRES -DPOSRES_LIG ; position restrain the protein and
> ligand
> ; Run parameters
> integrator = md ; leap-frog integrator
> nsteps = 50000 ; 2 * 50000 = 100 ps
> dt = 0.002 ; 2 fs
> ; Output control
> nstxout = 500 ; save coordinates every 1.0 ps
> nstvout = 500 ; save velocities every 1.0 ps
> nstenergy = 500 ; save energies every 1.0 ps
> nstlog = 500 ; update log file every 1.0 ps
> energygrps = Protein NPD
> ; Bond parameters
> continuation = no ; first dynamics run
> constraint_algorithm = lincs ; holonomic constraints
> constraints = hbonds ; all bonds (even heavy atom-H bonds)
> constrained
> lincs_iter = 1 ; accuracy of LINCS
> lincs_order = 4 ; also related to accuracy
> ; Neighborsearching
> cutoff-scheme = Verlet
> ns_type = grid ; search neighboring grid cells
> nstlist = 10 ; 20 fs, largely irrelevant with Verlet
> rcoulomb = 1.4 ; short-range electrostatic cutoff (in nm)
> rvdw = 1.4 ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype = PME ; Particle Mesh Ewald for long-range
> electrostatics
> pme_order = 4 ; cubic interpolation
> fourierspacing = 0.16 ; grid spacing for FFT
> ; Temperature coupling
> tcoupl = V-rescale ; modified Berendsen thermostat
> tc-grps = Protein_NPD Water_and_ions ; two coupling groups - more
> accurate
> tau_t = 0.1 0.1 ; time constant, in ps
> ref_t = 300 300 ; reference temperature, one
> for each group, in K
> ; Pressure coupling
> pcoupl = no ; no pressure coupling in NVT
> ; 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 = 300 ; temperature for Maxwell distribution
> gen_seed = -1 ; generate a random seed
>
> It generated no error in the log file. This is em.log
> Steepest Descents converged to Fmax < 100 in 10536 steps
> Potential Energy = -2.1753502e+06
> Maximum force = 8.4056786e+01 on atom 100002
> Norm of force = 4.0102091e+00
>
> For the NPT.mdp, since I have got error from the previous run, so I tried
> commenting out the neighborsearching part to take the default number and
> there is still error as the previous email.
>
This doesn't make sense to me. Please copy and paste error messages rather than
trying to remember or reinterpret them. Changing the nonbonded setup midway
through a run (or between stages) is fundamentally unsound. Verify that your
settings are right; IIRC a 1.4-nm cutoff is wrong for AMBER force fields.
-Justin
> Regarding the coordinates, I merged 2 coordinate files into 1 file as your
> tutorial about Protein-Ligand Simulation (
> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/complex/02_topology.html)
> and
> run a genconf to renumber the system, editconf to create the cubic box and
> genbox to solvate the system. Then I add 24 NA since the charge of the
> whole system is -24.
>
> I don't know which step is a possible cause to the error.
> Regards,
> Jennifer
>
>
> On Mon, Feb 2, 2015 at 1:57 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>
>>
>>
>> On 2/2/15 6:37 AM, Jennifer Vo wrote:
>>
>>> Dear Experts,
>>> I am facing the problem of simulation of protein - ligand complex using
>>> amber99sb force field.
>>> Since I created a topol.top for the system
>>>
>>> ; Include forcefield parameters
>>> #include "amber99sb.ff/forcefield.itp"
>>> [ atomtypes ]
>>> ;name bond_type mass charge ptype sigma epsilon
>>> Amb
>>> CA CA 0.00000 0.00000 A 3.39967e-01 3.59824e-01
>>> ;
>>> 1.91 0.0860
>>> H4 H4 0.00000 0.00000 A 2.51055e-01 6.27600e-02
>>> ;
>>> 1.41 0.0150
>>> HA HA 0.00000 0.00000 A 2.59964e-01 6.27600e-02
>>> ;
>>> 1.46 0.0150
>>> C C 0.00000 0.00000 A 3.39967e-01 3.59824e-01
>>> ;
>>> 1.91 0.0860
>>> O O 0.00000 0.00000 A 2.95992e-01 8.78640e-01
>>> ;
>>> 1.66 0.2100
>>> N N 0.00000 0.00000 A 3.25000e-01 7.11280e-01
>>> ;
>>> 1.82 0.1700
>>> H H 0.00000 0.00000 A 1.06908e-01 6.56888e-02
>>> ;
>>> 0.60 0.0157
>>> N* N* 0.00000 0.00000 A 3.25000e-01 7.11280e-01
>>> ;
>>> 1.82 0.1700
>>> CT CT 0.00000 0.00000 A 3.39967e-01 4.57730e-01
>>> ;
>>> 1.91 0.1094
>>> H2 H2 0.00000 0.00000 A 2.29317e-01 6.56888e-02
>>> ;
>>> 1.29 0.0157
>>> H1 H1 0.00000 0.00000 A 2.47135e-01 6.56888e-02
>>> ;
>>> 1.39 0.0157
>>> OH OH 0.00000 0.00000 A 3.06647e-01 8.80314e-01
>>> ;
>>> 1.72 0.2104
>>> HO HO 0.00000 0.00000 A 0.00000e+00 0.00000e+00
>>> ;
>>> 0.00 0.0000
>>> OS OS 0.00000 0.00000 A 3.00001e-01 7.11280e-01
>>> ;
>>> 1.68 0.1700
>>> P P 0.00000 0.00000 A 3.74177e-01 8.36800e-01
>>> ;
>>> 2.10 0.2000
>>> O2 O2 0.00000 0.00000 A 2.95992e-01 8.78640e-01
>>> ;
>>> 1.66 0.2100
>>> CK CK 0.00000 0.00000 A 3.39967e-01 3.59824e-01
>>> ;
>>> 1.91 0.0860
>>> H5 H5 0.00000 0.00000 A 2.42146e-01 6.27600e-02
>>> ;
>>> 1.36 0.0150
>>> NB NB 0.00000 0.00000 A 3.25000e-01 7.11280e-01
>>> ;
>>> 1.82 0.1700
>>> CB CB 0.00000 0.00000 A 3.39967e-01 3.59824e-01
>>> ;
>>> 1.91 0.0860
>>> N2 N2 0.00000 0.00000 A 3.25000e-01 7.11280e-01
>>> ;
>>> 1.82 0.1700
>>> NC NC 0.00000 0.00000 A 3.25000e-01 7.11280e-01
>>> ;
>>> 1.82 0.1700
>>> CQ CQ 0.00000 0.00000 A 3.39967e-01 3.59824e-01
>>> ;
>>> 1.91 0.0860
>>> ; Include chain topologies
>>> #include "A.itp"
>>> ; Include Position restraint file
>>> #ifdef POSRES
>>> #include "posre_A.itp"
>>> #endif
>>> ; Include chain topologies
>>> #include "B.itp"
>>> ; Include Position restraint file
>>> #ifdef POSRES
>>> #include "posre_B.itp"
>>> #endif
>>> ; Include custom ligand topologies
>>> #include "npd.itp"
>>> ; Ligand position restraints
>>> #ifdef POSRES_LIG
>>> #include "posre_npd.itp"
>>> #endif
>>>
>>> ; Include water topology
>>> #include "amber99sb.ff/tip3p.itp"
>>>
>>> #ifdef POSRES_WATER
>>> ; Position restraint for each water oxygen
>>> [ position_restraints ]
>>> ; i funct fcx fcy fcz
>>> 1 1 1000 1000 1000
>>> #endif
>>> ; Include generic ion topology
>>> #include "amber99sb.ff/ions.itp"
>>>
>>> [ system ]
>>> ; Name
>>> system in water
>>> [ molecules ]
>>> ; Compound #mols
>>> A 1
>>> B 1
>>> NPD 1
>>> SOL 42322
>>> NA 24
>>>
>>> The Protein and the Ligand were under energy minimization and a short MD
>>> separatedly to be sure there is no problem with the topology. But when I
>>> merge two topol.top from protein and from ligand, the error came at NPT
>>> step after successfully energy minimization and NVT
>>> "This is likely either a 1,4 interaction, or a listed interaction inside a
>>> smaller molecule you are decoupling during a free energy calculation.
>>> Since
>>> interactions at distances beyond the table cannot be computed, they are
>>> skipped until they are inside the table limit again. You will only see
>>> this
>>> message once, even if it occurs for several interactions.
>>> IMPORTANT: This should not happen in a stable simulation, so there is
>>> probably something wrong with your system. Only change the table-extension
>>> distance in the mdp file if you are really sure that is the reason."
>>>
>>> This is NPT.mdp
>>> title = system in water
>>> define = -DPOSRES -DPOSRES_LIG ; position restrain the protein
>>> ; Run parameters
>>> integrator = md ; leap-frog integrator
>>> nsteps = 500000 ; 2 * 500000 = 1000 ps, 1 ns
>>> dt = 0.002 ; 2 fs
>>> ; Output control
>>> nstxout = 5000 ; save coordinates every 1.0 ps
>>> nstvout = 5000 ; save velocities every 1.0 ps
>>> nstenergy = 500 ; save energies every 1.0 ps
>>> nstlog = 500 ; update log file every 1.0 ps
>>> energygrps = Protein NPD
>>> ; Bond parameters
>>> continuation = yes ; Restarting after NVT
>>> constraint_algorithm = lincs ; holonomic constraints
>>> constraints = all-bonds ; all bonds (even heavy atom-H bonds)
>>> constrained
>>> lincs_iter = 1 ; accuracy of LINCS
>>> lincs_order = 4 ; also related to accuracy
>>>
>>> ; Neighborsearching
>>> ;ns_type = grid ; search neighboring grid cels
>>> ;nstlist = 5 ; 10 fs
>>> ;rlist = 1.4 ; short-range neighborlist cutoff (in nm)
>>> ;rcoulomb = 1.4 ; short-range electrostatic cutoff (in
>>> nm)
>>> ;rvdw = 1.4 ; short-range van der Waals cutoff (in
>>> nm)
>>>
>>>
>> You have commented out the most critical elements of the input file. You
>> need to specify your nonbonded settings, otherwise you're relying solely on
>> defaults, which may or may not be right for the force field you're using.
>>
>>
>> ; Electrostatics
>>> coulombtype = PME ; Particle Mesh Ewald for
>>> long-range electrostatics
>>> pme_order = 4 ; cubic interpolation
>>> fourierspacing = 0.16 ; grid spacing for FFT
>>>
>>> ; Temperature coupling is on
>>> tcoupl = V-rescale ; modified Berendsen
>>> thermostat
>>> tc-grps = Protein_NPD Water_and_ions ; two coupling groups -
>>> more accurate
>>> tau_t = 0.1 0.1 ; time constant, in ps
>>> ref_t = 300 300 ; reference temperature, one for
>>> each group, in K
>>> ; Pressure coupling is on
>>> pcoupl = Berendsen ; Pressure coupling on in 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
>>> refcoord_scaling = com
>>> ; Periodic boundary conditions
>>> pbc = xyz ; 3-D PBC
>>> ; Dispersion correction
>>> DispCorr = EnerPres ; account for cut-off vdW scheme
>>> ; Velocity generation
>>> gen_vel = no ; Velocity generation is off
>>>
>>>
>> Given that you're not generating velocities, presumably you're picking up
>> from a previous run? What were its settings and ensemble? Was it stable?
>>
>> Also be sure to check the merged coordinates and make sure energy
>> minimization was sensible.
>>
>> -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
>>
>> ==================================================
>> --
>> Gromacs Users mailing list
>>
>> * Please search the archive at http://www.gromacs.org/
>> Support/Mailing_Lists/GMX-Users_List before posting!
>>
>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>> * For (un)subscribe requests visit
>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>> send a mail to gmx-users-request at gromacs.org.
>>
--
==================================================
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
==================================================
More information about the gromacs.org_gmx-users
mailing list