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

[Mark Abraham]

Hi,

That could only happen if your grompp -c input used a different atom
ordering from that implied by your [system] section of your .top. grompp
warns about this, but maybe you didn't notice... (Or the structure you
loaded into VMD is a mismatch to the trajectory, so its heuristics for
guessing where bonds actually are get double-crossed.)

Mark

On Sat, Jul 11, 2015 at 1:15 AM Nathan K Houtz <nhoutz at purdue.edu> wrote:

> Actually, I think I found the problem. When i looked at the pdb files in
> vmd the first time I missed it, but a colleague had a hunch and I found
> that he was right! I think that shake is constraining atoms to the wrong
> molecules. Here's a screenshot I took where you can see what I mean:
> http://imgur.com/1bgThDv
>
> I'm still not sure how to fix it though. The molecule numbering in my .gro
> file is correct, and there are only 4 atoms per molecule (or 3 plus the
> dummy). How can I tell shake to look at each molecule individually? Or have
> I made a mistake somewhere? Here are my most recent files:
>
> .gro file: (only part of it -- there is a character limit. this is the
> first and last ten molecules with the box dimensions at the bottom)
> http://textuploader.com/e898
>
> topology file:
> http://textuploader.com/e89g
>
> minimization:
> http://textuploader.com/e8rc
>
> .mdp file:
> http://textuploader.com/e896
>
> log file from a failed simulation:
> http://textuploader.com/e80y
>
> Thanks for your help! If you do still want the full .gro file to run the
> simulation, I can look for a different file host.
> Nathan
>
>
> ----- Original Message -----
> From: "Justin Lemkul" <jalemkul at vt.edu>
> To: gmx-users at gromacs.org
> Sent: Friday, July 10, 2015 5:06:34 PM
> Subject: Re: [gmx-users] trouble diagnosing a simulation that's blowing up
> (shake not converging --> segmentation fault)
>
>
>
> On 7/9/15 9:23 PM, Nathan K Houtz wrote:
> > Thanks for your explanations, Dr. Lemkul.
> >
> > I had already corrected a couple of the things you suggested. Gromacs
> won't actually let me run with Nose-Hoover and Parrinello-Rahman together
> (or at least, it gives a warning not to do that and stops). I'd like to run
> in NVT anyway, so I had set Pcoupl to 'no'. I had also increased the
> cutoffs and my nstlist is 20.
> >
> > I have now also changed the tcoupl to v-rescale, but unfortunately that
> alone didn't help. So I'm not sure exactly what I can inspect to find the
> source of my error. I know that temperature, pressure, and total energy are
> all warning signs of bad things if they misbehave, but they were all fairly
> constant throughout the simulation. (This is for the flexible-model
> simulation:) The starting energy for 1700 water molecules at 120K was
> -1.08744e+05, and it finished at -1.06047e+05 with no significant
> variations on any step. Pressure and temperature were also fine. I
> attempted to look at the results visually in vmd but I might have done
> something wrong because the .trr file has some error in it and vmd crashes.
> But the starting geometry after minimization looks fine: none of the
> molecules were moved out of their position in the crystal.
> >
> > For the constrained molecules, the energy stays about the same (about
> -1e+5) for the first 8 timesteps and then quickly blows up to +2e15 at step
> 12 before it obviously fails. I get shake warnings from step 0 though. In
> vmd, the first 8 steps look reasonable (checking the .pdb files that
> gromacs outputs when there are warnings) and all the molecules seem to hold
> their places in the crystal, but then at step 9 suddenly some of the
> molecules become misshapen with very long or very short bonds. Of course it
> goes downhill from there, but I can't figure out what's causing the
> problems since it's clearly happening from the very beginning (according to
> the shake warnings).
> >
> > What other things could I look at to troubleshoot my problem?
> >
>
> Not sure, but at least the failure happens fast.  Upload all your input
> files
> somewhere and I'll take a few minutes to see if I can spot anything.
>
> -Justin
>
> > Thanks for your help,
> > Nathan
> >
> > ----- Original Message -----
> > From: "Justin Lemkul" <jalemkul at vt.edu>
> > To: gmx-users at gromacs.org
> > Sent: Thursday, July 9, 2015 7:39:47 AM
> > Subject: Re: [gmx-users] trouble diagnosing a simulation that's blowing
> up (shake not converging --> segmentation fault)
> >
> >
> >
> > 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
>
> ==================================================
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>