Re: [AMBER-Developers] FW: [AMBER] Issue with running sander with SHAKE

From: Duke, Robert E Jr <>
Date: Wed, 7 Dec 2011 03:18:05 +0000

By the way, I am pretty sure pmemd just defines molecules completely differently, in order to get away from the ordering problem. What pmemd can't do for you however, is conceptualize things that are not somehow covalently bonded as one molecule. This seemed to me as close as correct as pmemd could get, given that you could not place much faith in pdb order. I understand the possibility of strongly ionically bonded exceptions...
Regards - Bob

From: Jason Swails []
Sent: Tuesday, December 06, 2011 7:13 PM
To: AMBER Developers Mailing List
Cc: Breuer, Marian
Subject: Re: [AMBER-Developers] FW: [AMBER] Issue with running sander with SHAKE

On Tue, Dec 6, 2011 at 8:43 PM, case <> wrote:

> On Tue, Dec 06, 2011, Breuer, Marian wrote:
> On the to-do list (increasing level of difficulty):
> 1. LEaP definitely should check for this, and certainly should not create
> an ATOMS_PER_MOLECULE data structure that is incorrect.

Agreed here. However, I know I have seen a tleap-created topology file
whose ATOMS_PER_MOLECULE and SOLVENT_POINTERS were so wrong that even pmemd
choked on it. The only thing I know about it is that it's an unusual
system to say the least. At the very least this shows there are some
corner cases that tleap doesn't handle correctly.

> At an even more
> basic minimum, Dave Cerutti's [or Jason's] prmtop validation/editing
> program
> could also flag this problem.]

My program does check for it, and also offers an option to recompute a
(completely correct) ATOMS_PER_MOLECULE and SOLVENT_POINTERS section,
choosing whether to set ions as part of the solute (tleap and sometimes
sleap default) or as part of the solvent (sleap default if water is added
first). However, if molecules (by topology) are not continuous atoms in
the prmtop, setting ATOMS_PER_MOLECULE and SOLVENT_POINTERS exits with an
error rather than combining everything into a single molecule.

2. The Users' manual should document this requirement, and error messages
> about its violation should point to such documentation.
> 3. The pmemd sanity check on ATOMS_PER_MOLECULE should be backported to
> sander.

I thought it was. If not, I'll backport it. I put that check in place
after someone gave me the topology file mentioned above.

> 4. (Longer term) better barostats use an "atomic" model for volume
> fluctuations, avoiding the rather arbitrary division into molecules.
> >
> > But an additional reason for my slow reply: my memory of the code is that
> > the splitting of atoms should depend on residues, not molecules. So I
> want
> > to see if I can understand how your original problem (with > 14 nodes)
> arose.
> As per usual, my memory was wrong. For constant pressure runs, splitting
> of atoms among processors depends on molecules, not residues [as comments
> in the set.f file (that I think I wrote!) clearly indicate.] If the
> information in the prmtop file is wrong, things can be badly messed up.

Very true. In the example I've mentioned to date, not every atom was
assigned to a molecule (since sum(ATOMS_PER_MOLECULE) < NATOM), which meant
that the molecule owner array (pmemd) was filled with garbage for all atoms
not assigned in ATOMS_PER_MOLECULE. Since this is a pointer array, the
result was a segfault. Sander resulted in weird behavior (no segfault,
just garbage results), and differently in serial and parallel.

All the best,

Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
AMBER-Developers mailing list
AMBER-Developers mailing list
Received on Tue Dec 06 2011 - 19:30:04 PST
Custom Search