Dear Users:
Summary:
I believe that there is some problem in the residue numbering routine that is
used by editconf and trjconv
that exists at least in 4.5.3 and 4.5.5, but not in 4.0.7, and that this is
somehow related to the ACE residue.
Details:
I have a system with two peptides, water, and ions. Residues are numbered
consecutively in the initial .gro file
(verified by running editconf -resnr 1). The problem is that when I use grompp
to generate a .tpr file and then
use editconf -f my.tpr -o out.gro I get a .gro file with residues that are not
listed consecutively (same thing if I use
trjconv to extract a .gro from a .xtc file with this .tpr).
I tested this in gromacs 4.5.5 and 4.5.3. The forcefield is either Amber99sb or
OPLS-AA/L. Each peptide is capped
(starts with ACE and ends with NME).
If I order them as protein1, protein2, water, ions in the .gro and .top then
the .gro that I extract from the .tpr
has residue numbers that run from 1-14 (first peptide), then again 1-15 (second
peptide), and continue as
16-1413 (water), 1414-1439 (ions).
If I reorder to put ions first, then I get residue numbers that run from 16-41
(ions), 1-14 (first peptide),
1-15 (second peptide), 42-1439 (water).
In each case, however, the output from gmxdump provides molblock fields that
are in the correct order.
I wonder if there could be something about the ACE residue in the force field
that is causing this
strange behaviour? I've used a single capped peptide before, but never more
than one and the ACE residue
has always come first in my previous work.
Here is the molblock information from gmxdump in the setup in which I put the
ions before the first peptide:
topology:
name="MD"
#atoms = 4652
molblock (0):
moltype = 0 "NA"
#molecules = 15
#atoms_mol = 1
molblock (1):
moltype = 1 "CL"
#molecules = 11
#atoms_mol = 1
molblock (2):
moltype = 2 "Protein_1"
#molecules = 1
#atoms_mol = 201
molblock (3):
moltype = 3 "Protein_2"
#molecules = 1
#atoms_mol = 231
molblock (4):
moltype = 4 "SOL"
#molecules = 1398
#atoms_mol = 3
#############
And here is the first 30 lines of the input .gro
RUN1
4652
1NA NA 1 1.111 1.195 1.100
2NA NA 2 0.903 0.086 0.133
...
25CL CL 25 2.684
2.864<tel:25%C2%A0%C2%A0%202.684%C2%A0%C2%A0%202.864> 3.234
26CL CL 26 0.233
0.729<tel:26%C2%A0%C2%A0%200.233%C2%A0%C2%A0%200.729> 2.132
27ACE CH3 27 3.073 2.376 1.686
27ACE HH31 28 3.163 2.414 1.665
...
###############################
And here is the first 30 lines of the output .gro after grompp and editconf -f
.tpr
RUN1
4652
16NA NA 1 1.111 1.195 1.100
17NA NA 2 0.903 0.086 0.133
...
40CL CL 25 2.684
2.864<tel:25%C2%A0%C2%A0%202.684%C2%A0%C2%A0%202.864> 3.234
41CL CL 26 0.233
0.729<tel:26%C2%A0%C2%A0%200.233%C2%A0%C2%A0%200.729> 2.132
1ACE CH3 27 3.073 2.376 1.686
1ACE HH31 28 3.163 2.414 1.665
...
########################
Finally, I tested this by creating an OPLS topology for either ALA or ACE-ALA.
I then put two molecules
of the selected species in a box and ran grompp; editconf -f .tpr -o .gro and
found that the atom numbers
started again back at zero only for the system containing ACE. This narrows it
down to show that problem is
not system specific, nor does it depend on the presence of the NME residue.
Further, I tested with version 4.0.7 (pdb2gmx, grompp, editconf) and found that
this problem did not exist.
Nevertheless, if I used version 4.5.5 of editconf to extract from the .tpr
generated by 4.0.7 then I see improper
residue numbering.
It seems to me that there is some problem in the residue numbering routine that
is used by editconf and trjconv
that exists at least in 4.5.3 and 4.5.5 but not in 4.0.7. This can be
problematic for later analysis if, for example, a .ndx file is constructed from
a .gro with non-sequential residue numbers. It also seems to be leading to
problems with our selections for genion.
If anybody has encountered this and solved it, I'd love to hear about it.
Thank you,
Chris.
--
gmx-users mailing list gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
