[gmx-users] g_energy inconsistent results

[Mark Abraham]

On 8/03/2011 9:44 PM, Ehud Schreiber wrote:
>
> Dear Gromacs users,
>
> I am working with version 4.5.3, using the opls-aa forcefield in an 
> implicit solvent, all-vs-all setting:
>
> pdb2gmx -ter -ff oplsaa -water none -f file.pdb
>
> I am energy-minimizing structures in 3 stages (steep, cg and l-bfgs). 
> The last stage is the following:
>
> grompp -f em3.mdp -p topol.top -c em2.gro -t em2.trr -o em3.tpr -po 
> em3.mdout.mdp
>
> mdrun -nice 0 -v -pd -deffnm em3
>
> g_energy -s em3.tpr -f em3.edr -o em3.potential_energy.xvg
>
> where the mdp file is:
>
> ;;;;;;;;;;;;;;;;;;; em3.mdp ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> integrator       = l-bfgs
>
> nsteps           = 50000
>
> implicit_solvent = GBSA
>
> gb_algorithm     = Still
>
> sa_algorithm     = Ace-approximation
>
> pbc              = no
>
> rgbradii         = 0
>
> ns_type          = simple
>
> nstlist          = 0
>
> rlist            = 0
>
> coulombtype      = cut-off
>
> rcoulomb         = 0
>
> vdwtype          = cut-off
>
> rvdw             = 0
>
> nstcalcenergy    = 1
>
> nstenergy        = 1000
>
> emtol            = 0
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> The last line in the em3.potential_energy.xvg file should give the 
> (potential) energy of the minimized structure em3.gro .
>
> I wish also to compute the potential energy of .gro files in general, 
> not necessarily obtained from a simulation. For that, I prepared a 
> .mdp file for a degenerate energy minimization, having 0 steps, 
> designed just to give the status of the file:
>

Zero-step EM does not calculate the initial energy because it is not 
useful for gradient-based energy minimization. I don't recall the 
details, but perhaps the first EM step is reported as step zero.

Instead, you should use zero-step MD (with unconstrained_start = yes), 
or (for multiple single points) mdrun -rerun.

You will not necessarily reproduce the g_energy energies with anything, 
because the energy is dependent on the state of the neighbour lists. If 
nstenergy is a multiple of nstlist, then those energies should be fairly 
reproducible.

I have updated the grompp source to provide a note to the user to warn 
against zero-step EM.

Mark

> ;;;;;;;;;;;;;;;;;;; status.mdp ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> integrator       = l-bfgs
>
> nsteps           = 0
>
> implicit_solvent = GBSA
>
> gb_algorithm     = Still
>
> sa_algorithm     = Ace-approximation
>
> pbc              = no
>
> rgbradii         = 0
>
> ns_type          = simple
>
> nstlist          = 0
>
> rlist            = 0
>
> coulombtype      = cut-off
>
> rcoulomb         = 0
>
> vdwtype          = cut-off
>
> rvdw             = 0
>
> nstcalcenergy    = 1
>
> nstenergy        = 1
>
> emtol            = 0
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> The only changes from the former .mdp file are in nsteps and nstenergy.
>
> However, when I run this potential energy status run on em3.gro itself,
>
> grompp -f status.mdp -p topol.top -c em3.gro -o status.tpr -po 
> status.mdout.mdp
>
> mdrun -nice 0 -v -pd -deffnm status
>
> g_energy -s status.tpr -f status.edr -o status.potential_energy.xvg
>
> and look at the (single) energy line in status.potential_energy.xvg I 
> find that the energy does not agree with the one obtained during 
> minimization (it's higher by some tens of kJ/mol).
>
> What am I doing wrong? How should one reliably find the energy of a 
> given .gro file?
>
> Moreover, when changing in status.mdp to integrator = steep, the 
> results also change dramatically -- why should the algorithm matter if 
> no steps are performed and the initial structure is explored?
>
> Thanks,
>
> Ehud.
>

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20110308/3789d11b/attachment.html>