[gmx-users] Energy Minimisation and Equilibration Problems

Nancy nancy5villa at gmail.com
Thu Aug 6 22:08:18 CEST 2009


I have attempted to perform energy minimisation from scratch again, and
these are the commands that I am using to do so:

$ .../topolbuild1_2_1/src/topolbuild -n ethanediol -dir
.../topolbuild1_2_1/dat/gromacs -ff gmx53a6

which generates the files:

ethanediol.gro
ethanediol.log
ethanediolMOL.mol2
ethanediol.top
ffethanediol.itp
posreethanediol.itp

in the "ethanediol.log" file I noticed the following lines with asterisks:

========================================
    Angles Force Field Results
Angle        Atoms        force    angle     method measured
   1   H12-   C1-  H11   ******   ******   ******    109.403
   2    C2-   C1-  H11   ******   ******   ******    109.456
   3    O1-   C1-  H11   ******   ******   ******    109.438
   4    C2-   C1-  H12   ******   ******   ******    109.537
   5    O1-   C1-  H12   ******   ******   ******    109.468
   6    O1-   C1-   C2    320.0   109.50        1    109.526
   7    C1-   C2-   O2    320.0   109.50        1    109.526
   8    C1-   C2-  H22   ******   ******   ******    109.453
   9    C1-   C2-  H21   ******   ******   ******    109.537
  10   H22-   C2-   O2   ******   ******   ******    109.439
  11   H21-   C2-   O2   ******   ******   ******    109.468
  12    C2-   O2-  HO2    450.0   109.50        1    106.864
  13   H21-   C2-  H22   ******   ******   ******    109.404
  14    C1-   O1-  HO1    450.0   109.50        1    106.864
========================================

I then enlarge the box and solvate the molecule:

$ editconf -f ethanediol.gro -o ethanediol_box.gro -box 3 3 3 -c

$ genbox -cp ethanediol_box.gro -cs spc216.gro -o ethanediol_solv.gro -p
ethanediol.top -box 3 3 3

There are 2679 atoms (891 waters and 1 ethanediol).  I then configure the
minimisation with grompp:

$ grompp -f minim.mdp -c ethanediol_solv.gro -p ethanediol.top -o em.tpr

using the following .mdp file:

========================================
define        = -DFLEXIBLE
integrator    = steep
emtol        = 10.0
emstep          = 0.1
nsteps        = 15000

nstlist        = 1
ns_type        = grid
rlist        = 1.0
coulombtype    = PME
rcoulomb_switch    = 1.0
rvdw_switch    = 1.3
pbc        = xyz

pme_order    = 4
constraints     = none

nstxout        = 1
nstvout        = 1
nstenergy    = 1
nstlog        = 1
nstcomm         = 1
Tcoupl          = no
Pcoupl          = no
gen_vel         = no
========================================

I then run the simulation for 15,000 steps:

$ mdrun -v -deffnm em

the last few lines of output are as follows:

========================================
Step=14991, Dmax= 5.3e-06 nm, Epot= -6.24538e+04 Fmax= 2.50258e+04, atom=
2397
Step=14992, Dmax= 6.4e-06 nm, Epot= -6.24539e+04 Fmax= 3.52982e+04, atom=
1171
Step=14993, Dmax= 7.6e-06 nm, Epot= -6.24540e+04 Fmax= 3.71731e+04, atom=
2397
Step=14994, Dmax= 9.2e-06 nm, Epot= -6.24540e+04 Fmax= 4.93409e+04, atom=
1171
Step=14995, Dmax= 1.1e-05 nm, Epot= -6.24541e+04 Fmax= 5.53109e+04, atom=
2397
Step=14997, Dmax= 6.6e-06 nm, Epot= -6.24546e+04 Fmax= 8.50832e+03, atom=
1171
Step=14998, Dmax= 7.9e-06 nm, Epot= -6.24547e+04 Fmax= 6.37691e+04, atom=
2397
Step=14999, Dmax= 9.5e-06 nm, Epot= -6.24553e+04 Fmax= 2.67488e+04, atom=
1171
Step=15000, Dmax= 1.1e-05 nm, Epot= -6.24548e+04 Fmax= 8.09451e+04, atom=
2397
writing lowest energy coordinates.

Steepest Descents did not converge to Fmax < 10 in 15001 steps.
Potential Energy  = -6.2455254e+04
Maximum force     =  2.6748828e+04 on atom 1171
Norm of force     =  9.5851715e+02
========================================

After viewing the trajectory with ngmx, I noticed that most of the motion of
the water molecules occurs within the first 50 picoseconds.

Please advise on how to further minimise this system.

Thank you.

Nancy




On Thu, Aug 6, 2009 at 2:23 PM, Bruce D. Ray <brucedray at yahoo.com> wrote:

> On Thursday, August 6, 2009 at 11:58:20 AM, Nancy <nancy5villa at gmail.com>
> wrote:
> >  The energies simply do not seem to come down any further within several
> thousand steps.
> >
> > I start with a .mol2 file which contains the structure of ethylene glycol
> (ethanediol).
> > These are the commands that I use to set up and run the minimisation:
> >
> > $ .../topolbuild1_2_1/src/topolbuild -n ethanediol -dir
> .../topolbuild1_2_1/dat/gromacs -ff gmx53a6
> >
> > the above command outputs the following files:
> >
> > ethanediol.gro
> > ethanediol.log
> > ethanediolMOL.mol2
> > ethanediol.top
> > ffethanediol.itp
> > posreethanediol.itp
>
> Those are the correct files to be generated for the topolbuild command
> given.
> Are the parameters blank in any of the lines in ethanediol.top?  Are there
> any lines in ethanediol.log that do not involve nonpolar hydrogens that
> have asterisks instead of parameter values?
>
> > I then proceed to enlarge the box and solvate the molecule:
> >
> > $ editconf -f ethanediol.gro -o ethanediol_box.gro -box 5 5 5
> >
> > $ genbox -cp ethanediol_box.gro -cs spc216.gro -o ethanediol_box.gro -p
> ethanediol.top -shell 1
>
> Is this a shell of solvent around the molecule rather than a box of solvent
> containing
> the molecule of interest?  How is that to work with vacuum outside the
> shell?
>
> > I then use grompp to configure the minimisation:
> >
> > $ grompp -f minim.mdp -c ethanediol_box.gro -p ethanediol.top -o em.tpr
> >
> > This is my .mdp file for minimisation:
> >
> > ===============================
> > define        = -DPOSRE
>
> Why position restraints when doing energy minimization?
>
> > integrator    = steep
> > emtol        = 10.0
> > emstep          = 0.01
> > nsteps        = 2000
>
> Why so few steps when the emtol is this small.  For this emtol, I would
> expect
> to need 10000 to 20000 steps.
>
> > nstlist        = 1
> > ns_type        = grid
> > rlist        = 1.0
> > coulombtype    = PME
> > rcoulomb    = 1.0
> > rvdw        = 1.3
> > pbc        = xyz
> >
> > pme_order    = 4
> > constraints     = none
> >
> >
> > nstxout        = 1
> > nstvout        = 1
> > nstenergy    = 1
> > nstlog        = 1
> > nstcomm         = 1
> > Tcoupl          = no
> > Pcoupl          = no
> > gen_vel         = no
> > ===============================
>
> I've been told that the last three lines are not needed for energy
> minimization.
>
> > $ mdrun -v -deffnm em
> >
> > the minimisation runs without error and these are the last few lines of
> output:
> >
> > ===========================================
> > Step= 1990, Dmax= 1.1e-05 nm, Epot= -1.76427e+04 Fmax= 1.25524e+04, atom=
> 304
> > Step= 1991, Dmax= 1.3e-05 nm, Epot= -1.76433e+04 Fmax= 3.13881e+04, atom=
> 484
> > Step= 1992, Dmax= 1.6e-05 nm, Epot= -1.76443e+04 Fmax= 2.11785e+04, atom=
> 484
> > Step= 1993, Dmax= 1.9e-05 nm, Epot= -1.76445e+04 Fmax= 4.22489e+04, atom=
> 484
> > Step= 1994, Dmax= 2.3e-05 nm, Epot= -1.76456e+04 Fmax= 3.33021e+04, atom=
> 484
> > Step= 1996, Dmax= 1.4e-05 nm, Epot= -1.76468e+04 Fmax= 1.26110e+04, atom=
> 304
> > Step= 1997, Dmax= 1.6e-05 nm, Epot= -1.76474e+04 Fmax= 3.81619e+04, atom=
> 484
> > Step= 1998, Dmax= 2.0e-05 nm, Epot= -1.76486e+04 Fmax= 2.68240e+04, atom=
> 484
> > Step= 2000, Dmax= 1.2e-05 nm, Epot= -1.76496e+04 Fmax= 1.26501e+04, atom=
> 304
> >
> > writing lowest energy coordinates.
> >
> > Steepest Descents did not converge to Fmax < 10 in 2001 steps.
> > Potential Energy  = -1.7649566e+04
> > Maximum force     =  1.2650057e+04 on atom 30
> > Norm of force     =  1.7544460e+03
> > ===========================================
> >
> > If I proceed to equilibration after doing the above, the water molecules
> simply fly apart
> > (although not immediately).  Additionally, it seems that there are no
> interactions between
> > the waters during equilibration.
>
> Really no point in proceeding further until the energy minimization issues
> are solved.
>
> > I have tried to run the minimisation for a larger number of steps, but it
> does not help.
> > I have also tried to delete individual water molecules from the structure
> files, but doing
> > so simply causes the minimisation to "fixate" on another two molecules.
> I am not sure
> > what values of the energies are reasonable for this system, and how to
> minimise it further.
> > Please advise.
>
> I've never used the -shell with genbox and don't know what effect that
> would have on
> energy minimization.
>
> --
> Bruce D. Ray, Ph.D.
> Associate Scientist, and Operations Director
> NMR Center
> IUPUI
> Physics Dept.
> 402 N. Blackford St.
> Indianapolis, IN 46202-3273
>
>
>
> _______________________________________________
> gmx-users mailing list    gmx-users at 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-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20090806/7dcf19d8/attachment.html>


More information about the gromacs.org_gmx-users mailing list