Nancy wrote:
Hello,
I have energy minimised my ethylene glycol system. I try to run
position restrained equilibration with constant volume and temperature
on it with the following .mdp file:
=============================================
title = Ethanediol equilibration 1
define = -DPOSRES
; Run parameters
integrator = md
nsteps = 5000
dt = 0.0002
5000 steps * 0.0002 ps/step = 1 ps. I doubt that length of time is sufficient
to equilibrate even the simplest system.
<snip>
rlist = 1.0
rcoulomb = 1.0
rvdw = 1.3
If you're still using Gromos96, I'd suggest you read the original derivation of
the force field to get the correct values for these parameters, most notably rvdw.
<snip>
I ran the equilibration with the command:
$ grompp -f nvt.mdp -c em.gro -p ethanediol.top -o nvt.tpr
where "em.gro" is the output coordinates of the minimisation, and
"ethanediol.top" is the topology of the system. I then run the
equilibration with the following command:
You haven't told us much about the energy minimization. What Fmax and Epot did
your system converge to? Do you believe those values to be acceptable?
$ mdrun -v -deffnm nvt
The system was "blowing up", and mdrun would not complete the
equilbration. I tried decreasing dt, and although mdrun finishes
without error, the output trajectory (viewed through ngmx) shows that
the system has "blown up". Additionally, the structure of the
ethanediol solute seems to be incorrect throughout the equilbration.
Have you corrected for periodicity effects? The trajectory will look scrambled
if any molecules cross periodic boundaries.
Beyond this, I doubt anyone can provide any magical insights into what's going
wrong without a better description of what "incorrect" means in the context of
the structure.
-Justin
I am unsure what is going wrong, please help.
Thank you.
Nancy
On Wed, Aug 5, 2009 at 9:13 AM, Bruce D. Ray <bruced...@yahoo.com
<mailto:bruced...@yahoo.com>> wrote:
On Tuesday, August 4, 2009, at 6:30:57 PM, Nancy
<nancy5vi...@gmail.com <mailto:nancy5vi...@gmail.com>> wrote:
> I am trying to simulate a simple system of ethylene glycol
(ethane-1,2-diol) solvated in a water box.
> I converted the pdb file to a mol2 file and used topolbuild 1.2.1
to generate topology files:
>
> $ topolbuild -n ethanediol -dir .../topolbuild1_2_1/dat/gromacs
-ff gmx53a6
I prefer to use an absolute reference to the parameters directory, but
I also have my molecules files in directories outside of the topolbuild
directory.
> I the used editconf to enlarge the box. I then solvated the
molecule using the following command:
I presume that you checked the topology produced to be sure that
parameters
were found for all atoms, bonds, angles, dihedrals, and impropers in
the molecule
before proceeding. This is always a necessary step with any
automatic topology
generation program because the force field chosen might not have
parameters
for everything in your molecule. (In this case, I believe that
ethylene glycol
should not present any problems, but always checking the topology
generated
is a good habit to develop.)
> $ genbox -cp ethanediol_box.gro -cs spc216.gro -o
ethanediol_solv.gro -box 3 3 3 -p ethanediol.top
>
> where "ethanediol_box.gro" is the structure of the molecule. I
used grompp to generate the mdrun input file:
>
> $ grompp -f grompp.mdp -c ethanediol_solv.gro -p ethanediol.top
-o ethanediol.tpr
>
> grompp.mdp is the following:
>
> ========================================
> title = Ethanediol
> cpp = /lib/cpp
> include = -I../top
> define =
> integrator = md
> dt = 0.002
> nsteps = 50000
> nstxout = 1
> nstvout = 5
> nstlog = 5
> nstenergy = 10
> nstxtcout = 1
> xtc_grps = EDO SOL
> energygrps = EDO SOL
> nstlist = 10
> ns_type = grid
> rlist = 0.8
rlist seems short to me. Probably ought to use at least 1.0
> coulombtype = cut-off
Its best to use PME with topolbuild generated topologies.
> rcoulomb = 1.4
> rvdw = 0.8
rvdw is too small. Probably ought to be 1.2 to 1.4
> tcoupl = Berendsen
> tc-grps = EDO SOL
> tau_t = 0.1 0.1
> ref_t = 300 300
> Pcoupl = Berendsen
> tau_p = 1.0
> compressibility = 4.5e-5
> ref_p = 1.0
> gen_vel = yes
> gen_temp = 300
> gen_seed = 173529
> constraints = all-bonds
> ========================================
>
> where "EDO" refers to ethylene glycol. grompp creates several notes:
>
> NOTE 1 [file grompp.mdp, line unknown]:
> The Berendsen thermostat does not generate the correct kinetic
energy
> distribution. You might want to consider using the V-rescale
thermostat.
Something to think about.
> NOTE 2 [file grompp.mdp, line unknown]:
> You are using a plain Coulomb cut-off, which might produce
artifacts.
> You might want to consider using PME electrostatics.
>
> NOTE 3 [file grompp.mdp, line unknown]:
> This run will generate roughly 2500 Mb of data
That is a large amount of data. Perhaps make nstxout larger.
Also, this looks like a production simulation. I presume you did
an energy minimization step and a position restrained equilibration
before this, but you did not show us anything about these. Did they
have any errors?
> I then run the simulation:
>
> $ mdrun -s ethanediol.tpr -o traj.trr -x traj.xtc -v
>
> However mdrun produces several errors and warnings:
>
> t = 0.036 ps: Water molecule starting at atom 2659 can not be
settled.
> Check for bad contacts and/or reduce the timestep.
> Wrote pdb files with previous and current coordinates
>
> t = 0.038 ps: Water molecule starting at atom 2659 can not be
settled.
> Check for bad contacts and/or reduce the timestep.
> Wrote pdb files with previous and current coordinates
>
> t = 0.040 ps: Water molecule starting at atom 2659 can not be
settled.
> Check for bad contacts and/or reduce the timestep.
> Wrote pdb files with previous and current coordinates
>
> Step 21 Warning: pressure scaling more than 1%, mu: 1.0372
1.0372 1.0372
>
> Step 29, time 0.058 (ps) LINCS WARNING
> relative constraint deviation after LINCS:
> rms 2.714709, max 3.541824 (between atoms 2 and 3)
> bonds that rotated more than 30 degrees:
> atom 1 atom 2 angle previous, current, constraint length
> 1 4 88.8 0.6606 0.6855 0.1520
> 2 3 90.6 0.4305 0.4542 0.1000
> 1 2 89.2 0.6123 0.4690 0.1435
> 5 6 93.4 0.4193 0.1741 0.1000
> 4 5 89.2 0.6156 0.5033 0.1435
> ===================================
>
It is blowing up. Did you see anything in the energy minimization
or the
position restrained equilibration steps that seemed unusual?
>
> I am not sure where the errors are occuring (I am also wondering
if the absence of
> explicit non-polar hydrogens in the .gro files are relevant).
Please adivse.
Gromacs force fields are united atom force fields. As such,
explicit non-polar hydrogens
are supposed to be removed.
I hope these few remarks help.
Sincerely,
--
Bruce D. Ray, Ph.D.
Associate Scientist
IUPUI
Physics Dept.
402 N. Blackford St.
Indianapolis, IN 46202-3273
_______________________________________________
gmx-users mailing list gmx-users@gromacs.org
<mailto:gmx-users@gromacs.org>
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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
<mailto:gmx-users-requ...@gromacs.org>.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
------------------------------------------------------------------------
_______________________________________________
gmx-users mailing list gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
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
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
