# [gmx-users] Energy Minimisation and Equilibration Problems

Nancy nancy5villa at gmail.com
Fri Aug 7 01:59:43 CEST 2009

```I ran the minimisation, and mdrun gave the following last few lines of
output:

==================================
Step=19990, Dmax= 7.8e-06 nm, Epot= -7.14922e+04 Fmax= 3.19585e+04, atom=
2395
Step=19992, Dmax= 4.7e-06 nm, Epot= -7.14924e+04 Fmax= 8.72456e+03, atom=
2395
Step=19993, Dmax= 5.6e-06 nm, Epot= -7.14927e+04 Fmax= 3.70693e+04, atom=
2395
Step=19994, Dmax= 6.8e-06 nm, Epot= -7.14929e+04 Fmax= 2.10513e+04, atom=
2395
Step=19995, Dmax= 8.1e-06 nm, Epot= -7.14929e+04 Fmax= 4.71234e+04, atom=
2395
Step=19996, Dmax= 9.7e-06 nm, Epot= -7.14931e+04 Fmax= 3.68827e+04, atom=
2395
Step=19998, Dmax= 5.8e-06 nm, Epot= -7.14933e+04 Fmax= 1.35461e+04, atom=
2395
Step=19999, Dmax= 7.0e-06 nm, Epot= -7.14933e+04 Fmax= 4.77130e+04, atom=
2395
Step=20000, Dmax= 8.4e-06 nm, Epot= -7.14936e+04 Fmax= 2.41860e+04, atom=
2395

writing lowest energy coordinates.

Steepest Descents did not converge to Fmax < 10 in 20001 steps.
Potential Energy  = -7.1493609e+04
Maximum force     =  2.4185994e+04 on atom 2395
Norm of force     =  8.1511212e+02
==================================

As you can be seen, the forces still do not converge to Fmax < 10, even
after 20,000 steps.

Does anyone know what the problem might be?

Thanks,

Nancy

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

> On Thursday, August 6, 2009 at 4:08:18 PM, Nancy <nancy5villa at gmail.com>
> wrote:
>
> > 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
> > ========================================
>
> These are fine.  Only the nonpolar hydrogens lack parameters and they are
> removed from
> consideration in the final topology because this is a united atoms model
> force field.
>
> > 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
>
> Just out of curiosity, why are you defining flexible waters at this point?
>
>
> > 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.
> >
>
>
> I downloaded NSC93876 ethylene glycol as a Sybyl (.sy2) file from NCI,
> processed it with
> dos2unix, added the missing@<TRIPOS> last line, and named it
> ethylene_glycol.mol2
> Charge information is absent on this file as downloaded, and I did not do
> anything to
> account for charges.  This may make my results less than useful.  I then
> ran the following commands:
>
>      /usr/bdray/bin/topolbuild1_3 -dir /usr/bdray/tables -ff gmx53a6 -n
> glycol/ethylene_glycol
>      editconf -f ethylene_glycol.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
>      grompp -f em.mdp -c ethanediol_solv.gro -p ethanediol.top -o em.tpr
>
> File em.mdp includes:
>      constraints         =  none
>      integrator          =  steep
>      nsteps              =  20000
>      ;
>      ;    Energy minimizing stuff
>      ;
>      emtol               =  10.0
>      emstep              =  0.01
>
>      nstlist             =  1
>      coulombtype         =  PME
>      nstcomm             =  1
>      ns_type             =  grid
>      rlist               =  1
>      rcoulomb            =  1.0
>      rvdw                =  1.3
>      nstxout             =  1
>      pbc                 =  xyz
>      pme_order           =  4
>
> I deliberately did not use the flexible waters model with this energy
> minimization
> I then ran
>
>      mdrun -v -deffnm em
>
> I got the result:
>      Stepsize too small, or no change in energy.
>      Converged to machine precision,
>      but not to the requested precision Fmax < 10
>
>      Double precision normally gives you higher accuracy.
>      You might need to increase your constraint accuracy, or turn
>      off constraints alltogether (set constraints = none in mdp file)
>
>      Steepest Descents converged to machine precision in 6291 steps,
>      but did not reach the requested Fmax < 10.
>      Potential Energy  = -4.9487391e+04
>      Maximum force     =  9.0675652e+01 on atom 5
>      Norm of force     =  2.5599347e+02
>
> I note that atom 5 is part of the ethylene glycol molecule.   Furthermore,
> with very few
> exceptions, most of the activity in the energy minimization was devoted to
> the ethylene
> glycol molecule.  I presume that I would reach the desired Fmax if I were
> to use double
> precision.  However, an Fmax < 90.7 would seem to me to be usable.
>
> Perhaps the problem stems from the use of flexible waters in the energy
> minimization?
>
>
> --
> Bruce D. Ray, Ph.D.
> Associate Scientist
> 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.