[gmx-users] Energy Minimization with Lipid Bilayer

[Mark Abraham]

Justin Lemkul wrote:
>>>Steepest Descents converged to machine precision in 15 steps,
>>>but did not reach the requested Fmax < 1000.
>>>Potential Energy  =  5.9483964e+18
>>>Maximum force     =            inf on atom 6710
>>>Norm of force     =            nan
>>
>>How did you add these new waters? Should your system be periodic? Does
>>visualising the start, progress and end of the EM trajectory shed any
>>light on what is going wrong?
> 
> 
> The waters were added using a modification of the make_hole script that is
> widely used.  A member of our lab tweaked it a few years ago, and we apply it
> to a number of systems.  I created a box of water large enough to solvate my
> protein, and inserted the bilayer/crystal water system into it, with minimal
> overlap, but enough such that a few waters were deleted from the original box.

This sounds like you are adding a slab onto one end of the box, rather 
than replacing part of the box... did you adjust the box dimensions 
accordingly? (They're on the bottom line of a .gro file.) If so, gromacs 
might be getting confused about the positions of atoms outside the 
box... but this is me clutching at straws.

> Another note is that I have tried the minimization a couple times (with
> insignificant changes to my .mdp file, as it turns out), and the infinite force
> keeps changing location.  Is that maybe a clue as to what is going on?

Doesn't suggest anything to me.

> The system does not need to be periodic.  When I took a look at velocities and
> forces with g_traj (from the default output traj.trr), there were no velocities
> found (probably due to the fact that I did not set gen_vel = yes), and the
> forces were all 'nan.'  Is there another way to look at my trajectory output
> that might be of more use?  (Again, still new to Gromacs and very
> overwhelmed/impressed by the myriad function it can perform.)

Reality check - energy minimization does not involve velocities! See the 
introductory chapter(s) of the GROMACS manual if you don't understand why.

Use something like VMD (or the GROMACS built-in ngmx) to view the 
positions over the "trajectory" (EM provides a series of coordinate sets 
according to nstxout, not a true trajectory) - either the first step or 
final coordinates are probably enough to point out what is wrong. 
Viewing forces or velocities is a great way to start suffering from 
epilepsy.

> Thanks in advance for any more information you may be able to provide.

No problem - and thanks for the useful description.

Mark