[gmx-users] RE: g_energy inconsistent results

[Per Larsson]

Hi!

When you are computing your zero-step energies, do you then start from the gro-file that you got from em? If so, maybe the energies changes because gro-files have a fixed precision format (with three decimals), while the em-calculations are done using either full single or double precision.

Your second issue is almost certainly related to rounding errors. The all-vs-all and the cutoff-code (even with enormous cut-offs) compute interactions in completely order. I would not worry about a 0.005 nm difference in RMSD.

/Per

10 mar 2011 kl. 14.27 skrev Ehud Schreiber:

> Dear Mark Abraham (and anyone else interested),
> 
> I have implemented your suggestion, changing in the status.mdp file to
> integrator = md and adding continuation = yes (this is the new name of
> the unconstrained_start parameter); however, this did not have any
> effect - the energy remained as before, different from the one obtained
> during the minimization.
> 
> I have then encountered another (perhaps related) issue.
> I thought that the problem may arise from the combination of the
> generalized Born and all-vs-all settings.
> I have therefore changed in the minimization to rgbradii = rlist =
> rcoulomb = rvdw = 100 (from = 0).
> As my system is far smaller than 100 nm, I expected these parameters to
> provide also an all-vs-all setting (even if non-optimized).
> Nevertheless, the resulting structure differed from the previous
> minimized one (rmsd ~ 0.005 nm, delta energy ~ a few kJ/mol). Can this
> arise from rounding effects only or does this signify a problem? I'm not
> sure what rgbradii does, but the manual states that it must equal rlist.
> In addition, changing also status.mdp to have all radii = 100 and
> running on the new optimized output again does not yield the same energy
> as during the optimization. 
> 
> Does all of this make any sense to you?
> 
> Thank,
> Ehud.
> 
> 
>> Message: 6
>> Date: Tue, 08 Mar 2011 23:37:12 +1100
>> From: Mark Abraham <Mark.Abraham at anu.edu.au>
>> Subject: Re: [gmx-users] g_energy inconsistent results
>> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>> Message-ID: <4D7622F8.7050407 at anu.edu.au>
>> Content-Type: text/plain; charset="iso-8859-1"
>> 
>> 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.
>>> 
> --
> 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/Support/Mailing_Lists/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/Support/Mailing_Lists