Re: [AMBER] Issue with running sander with SHAKE

From: Breuer, Marian <marian.breuer.10.ucl.ac.uk>
Date: Sat, 5 Nov 2011 22:52:12 +0000

Dear Dr. Case, After some tests with my input files it seems that leap writes out incorrect ATOMS_PER_MOLECULE sections when different (i.e. TER-separated) molecules are connected by additional bonds using the "bond"-command. The hemes in my system are covalently bound to the protein via histidine ligands and cysteine linkages, and the necessary bonds I have to add myself upon loading the pdb file. Leaving out these bond commands leads to a correct ATOMS_PER_MOLECULE section. However, the bonds connect my ten hemes with two of the three protein chains (and these with each other), and the first number in the incorrect ATOMS_PER_MOLECULE section (8265) is the combined number of atoms in these twelve separate species. Following your suggestion, I've created a simple test system that reproduces this error - a sequence of four residues incl. two cysteines (CYX) in two separate chains. Simple loading the pdb file, setting the box and writing out a topology file using saveamberparmpol leads to a correct ATOMS_PER_MOLECULE section; however, if I manually bind the two cysteines before writing out the topology file the two subsets of atoms are added up again in the ATOMS_PER_MOLECULE section. By the way, without setting the box in leap using "set [...] box { ... }", no ATOMS_PER_MOLECULE section is written at all. Is it supposed to be like that? Is a different decomposition scheme used in cluster conditions that does not require the ATOMS_PER_MOLECULE section? I've attached the test input pdb; I used the following commands in xleap to create a topology file from it: tmp=loadpdb [Path]/CYX_only.pdb set tmp box { 110 100 110 } bond tmp.2.SG tmp.4.SG saveamberparmpol tmp [Path]/[Topologyfilename] [crdfilename] All the best, Marian ________________________________________ From: David A Case [case.biomaps.rutgers.edu] Sent: 04 November 2011 17:39 To: AMBER Mailing List Subject: Re: [AMBER] Issue with running sander with SHAKE On Fri, Nov 04, 2011, Breuer, Marian wrote: > > Thank you so much for this information. I used indeed an input pdb that > I had written out from xleap before; but I checked it and it did contain > all necessary TER specifiers. The hemes (incl. "HER") had already been > at the beginning above the protein atoms in the original pdb, and the > number of protein atoms written out by xleap (8365) is not even equal to > the actual number of protein atoms (8363; maybe xleap somehow counted > two TER lines in between the protein atom lines that terminate two > subchains of my protein). I thus would give parmed a try. It would help if you could post the input pdb file plus the exact commands you used to make the prmtop file. Or maybe the simplest system that shows the error (say with just one HEM group?)--it doesn't have to lead to a usable system, just an example that leads to bad pointers in the prmtop. [My best guess: there used to be some pieces of code in LEaP that ordered molecules by size (greatest number of atoms first, then go in descending order). I had thought that such behavior had all been removed, but something might be left. Since "real" PDB files almost always have heme groups and other ligands after the protein, instead of before it, there might be a bug lurking related to this. If so, putting the HEM groups after the protein might make a difference.] But if you can provide any example that doesn't work correctly, that is a great help in fixing the problem. Using parmed is only a work-around...we should try to track down the problem at its source. ...thx...dac _______________________________________________ AMBER mailing list AMBER.ambermd.org http://lists.ambermd.org/mailman/listinfo/amber

_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber

Received on Sat Nov 05 2011 - 16:00:02 PDT
Custom Search