[gmx-users] trouble diagnosing a simulation that's blowing up (shake not converging --> segmentation fault)

[Justin Lemkul]

On 7/8/15 8:06 PM, Nathan K Houtz wrote:
> Hello,
>
> I deleted the email and can't respond to my last reply directly - sorry! I
> got this response from Mark Abraham:
>
> Hi, Try doing some EM and initial equilibration with no constraints at all,
> perhaps? Mark
>
> I tried commenting out the shake commands, and got a short (5000 step)
> simulation to run just fine without blowing up. Before, I would get shake
> warnings from the first few steps and a segmentation fault around step 13 or
> 14. I would like to be able to simulate with rigid molecules, though. Why
> would the simulation work with flexible molecules but not rigid ones?
>

Flexible water allows weird geometry, which is probably coming up and causing 
your constraint algorithm to fail.  I'd inspect the outcome carefully.  Just 
because it runs doesn't mean it's right.

> Also, in the example .mdp file for tip4p water, there is the (outdated)
> option, 'unconstrained-start', which is now 'continuation'. I got errors when
> trying to make the input .tpr file when I attempted to set that option to
> 'yes'. The warning said it was because I want Gromacs to generate velocities
> to start the simulation, which is incompatible with that command. Is there
> another way I can try to start the simulation unconstrained? Or would you
> suggest another idea to fix my shake warnings?
>

What this setting says is "have the constraints already been solved 
(continuation = yes) or should mdrun constrain the starting configuration itself 
(continuation = no)."

The .mdp file has a number of weird settings.  I would never use Nose-Hoover and 
Parrinello-Rahman when generating velocities; that's likely to be very unstable. 
  See if a more forgiving thermostat and/or barostat resolves the issue, e.g. 
tcoupl = v-rescale and pcoupl = Berendsen.

Also note that the cutoffs are very short (probably because they want to use a 
small box, but note that it *does* affect the physics) and nstlist = 1 is 
totally unnecessary.  It doesn't hurt your physics, but it's a major waste of 
performance.

-Justin

-- 
==================================================

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

==================================================