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

From: Jason Swails <>
Date: Tue, 6 Dec 2011 22:13:14 -0500

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
Received on Tue Dec 06 2011 - 19:30:03 PST
Custom Search