<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body text="#000000" bgcolor="#ffffff">
On 8/03/2011 9:44 PM, Ehud Schreiber wrote:
<blockquote
cite="mid:64E0EFC300379B4EB6A055D519417854010340D0@cmail"
type="cite">
<meta http-equiv="Content-Type" content="text/html;
charset=ISO-8859-1">
<meta name="Generator" content="Microsoft Word 12 (filtered
medium)">
<style><!--
/* Font Definitions */
@font-face
        {font-family:Calibri;
        panose-1:2 15 5 2 2 2 4 3 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
        {margin:0in;
        margin-bottom:.0001pt;
        font-size:11.0pt;
        font-family:"Calibri","sans-serif";}
a:link, span.MsoHyperlink
        {mso-style-priority:99;
        color:blue;
        text-decoration:underline;}
a:visited, span.MsoHyperlinkFollowed
        {mso-style-priority:99;
        color:purple;
        text-decoration:underline;}
span.EmailStyle17
        {mso-style-type:personal-compose;
        font-family:"Calibri","sans-serif";
        color:windowtext;}
.MsoChpDefault
        {mso-style-type:export-only;}
@page WordSection1
        {size:8.5in 11.0in;
        margin:1.0in 1.25in 1.0in 1.25in;}
div.WordSection1
        {page:WordSection1;}
--></style><!--[if gte mso 9]><xml>
<o:shapedefaults v:ext="edit" spidmax="1026" />
</xml><![endif]--><!--[if gte mso 9]><xml>
<o:shapelayout v:ext="edit">
<o:idmap v:ext="edit" data="1" />
</o:shapelayout></xml><![endif]-->
<div class="WordSection1">
<p class="MsoNormal">Dear Gromacs users,<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">I am working with version 4.5.3, using the
opls-aa forcefield in an implicit solvent, all-vs-all setting:<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">pdb2gmx -ter -ff oplsaa -water none -f
file.pdb<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">I am energy-minimizing structures in 3
stages (steep, cg and l-bfgs). The last stage is the
following:<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">grompp -f em3.mdp -p topol.top -c em2.gro
-t em2.trr -o em3.tpr -po em3.mdout.mdp<o:p></o:p></p>
<p class="MsoNormal">mdrun -nice 0 -v -pd -deffnm em3<o:p></o:p></p>
<p class="MsoNormal">g_energy -s em3.tpr -f em3.edr -o
em3.potential_energy.xvg <o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">where the mdp file is:<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">;;;;;;;;;;;;;;;;;;; em3.mdp
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;<o:p></o:p></p>
<p class="MsoNormal">integrator = l-bfgs<o:p></o:p></p>
<p class="MsoNormal">nsteps = 50000<o:p></o:p></p>
<p class="MsoNormal">implicit_solvent = GBSA<o:p></o:p></p>
<p class="MsoNormal">gb_algorithm = Still <o:p></o:p></p>
<p class="MsoNormal">sa_algorithm = Ace-approximation <o:p></o:p></p>
<p class="MsoNormal">pbc = no<o:p></o:p></p>
<p class="MsoNormal">rgbradii = 0 <o:p></o:p></p>
<p class="MsoNormal">ns_type = simple<o:p></o:p></p>
<p class="MsoNormal">nstlist = 0<o:p></o:p></p>
<p class="MsoNormal">rlist = 0<o:p></o:p></p>
<p class="MsoNormal">coulombtype = cut-off<o:p></o:p></p>
<p class="MsoNormal">rcoulomb = 0<o:p></o:p></p>
<p class="MsoNormal">vdwtype = cut-off<o:p></o:p></p>
<p class="MsoNormal">rvdw = 0<o:p></o:p></p>
<p class="MsoNormal">nstcalcenergy = 1<o:p></o:p></p>
<p class="MsoNormal">nstenergy = 1000<o:p></o:p></p>
<p class="MsoNormal">emtol = 0 <o:p></o:p></p>
<p class="MsoNormal">;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">The last line in the
em3.potential_energy.xvg file should give the (potential)
energy of the minimized structure em3.gro .<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">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:</p>
</div>
</blockquote>
<br>
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.<br>
<br>
Instead, you should use zero-step MD (with unconstrained_start =
yes), or (for multiple single points) mdrun -rerun.<br>
<br>
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.<br>
<br>
I have updated the grompp source to provide a note to the user to
warn against zero-step EM.<br>
<br>
Mark<br>
<br>
<blockquote
cite="mid:64E0EFC300379B4EB6A055D519417854010340D0@cmail"
type="cite">
<div class="WordSection1">
<p class="MsoNormal"><o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">;;;;;;;;;;;;;;;;;;; status.mdp
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;<o:p></o:p></p>
<p class="MsoNormal">integrator = l-bfgs<o:p></o:p></p>
<p class="MsoNormal">nsteps = 0<o:p></o:p></p>
<p class="MsoNormal">implicit_solvent = GBSA<o:p></o:p></p>
<p class="MsoNormal">gb_algorithm = Still<o:p></o:p></p>
<p class="MsoNormal">sa_algorithm = Ace-approximation<o:p></o:p></p>
<p class="MsoNormal">pbc = no<o:p></o:p></p>
<p class="MsoNormal">rgbradii = 0<o:p></o:p></p>
<p class="MsoNormal">ns_type = simple<o:p></o:p></p>
<p class="MsoNormal">nstlist = 0<o:p></o:p></p>
<p class="MsoNormal">rlist = 0<o:p></o:p></p>
<p class="MsoNormal">coulombtype = cut-off<o:p></o:p></p>
<p class="MsoNormal">rcoulomb = 0<o:p></o:p></p>
<p class="MsoNormal">vdwtype = cut-off<o:p></o:p></p>
<p class="MsoNormal">rvdw = 0<o:p></o:p></p>
<p class="MsoNormal">nstcalcenergy = 1<o:p></o:p></p>
<p class="MsoNormal">nstenergy = 1<o:p></o:p></p>
<p class="MsoNormal">emtol = 0<o:p></o:p></p>
<p class="MsoNormal">;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">The only changes from the former .mdp file
are in nsteps and nstenergy.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">However, when I run this potential energy
status run on em3.gro itself,<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">grompp -f status.mdp -p topol.top -c
em3.gro -o status.tpr -po status.mdout.mdp<o:p></o:p></p>
<p class="MsoNormal">mdrun -nice 0 -v -pd -deffnm status<o:p></o:p></p>
<p class="MsoNormal">g_energy -s status.tpr -f status.edr -o
status.potential_energy.xvg<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">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).<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">What am I doing wrong? How should one
reliably find the energy of a given .gro file?<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">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?<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Thanks,<o:p></o:p></p>
<p class="MsoNormal">Ehud.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
</blockquote>
<br>
</body>
</html>