[gmx-users] dialanine using charmm27 ff
Sandeep Somani
ssomani at gmail.com
Sat Sep 10 02:23:40 CEST 2011
Hi Mark
It worked ! pdb2gmx is now properly processing the pdb, and moreover the
potential energy from gmx is within 1.27 kJ/mol of that from Charmm. More on
energy comparison shortly.
However, further tinkering of the rtp and pdb file was required to generate
the topology:
1.
The original rtp entry which had HT1/2/3 for the methyl hydrogens was giving
the error
"
Fatal error:
Atom H1 in residue ACE 1 was not found in rtp entry ACE with 6 atoms
while sorting atoms.
"
even though neither rtp nor pdb had H1 atom. I do not understand why it was
insisting on H1/2/3. It worked upon changing HT1/2/3 to H1/2/3.
2.
I changed bond entry 'N CAT' -> 'N CT' in rtp entry of CT3. Guess CAT was
a typo.
3.
The improper dihedral in CT3
' -C CT N -O ; commented to match with charmm psf '
was absent from the charmm psf. So I commented it.
fyi - final pdb and rtp entries are listed below.
Now, regarding the energies, I compared the different bonded and non-bonded
energy components from charmm with those from gmx.
Bonds, dihedral and coulomb energies matched to within 0.01 kJ/mol which was
great!
The discrepancy is in angles (U-B) where gmx value is higher by 1.6652
kJ/mol and LJ where gmx is lower by 0.459 kJ/mol.
Any idea on how to fix the angles (not too bothered by LJ) ?
I checked that the total number of angles in charmm psf and gmx top are the
same (=54), but havn't yet compared individual entries.
btw, would you (or someone else here) have a sample charmm (inp, crd, psf)
and gmx (mdp, gro, top) files where the two energies show a better match? It
would help me with this troubleshooting.
Thanks for your help
Sandeep
================final rtp entries============================
[ ACE ]
[ atoms ]
CT3 CT3 -0.270 0
H1 HA 0.090 1
H2 HA 0.090 2
H3 HA 0.090 3
C C 0.510 4
O O -0.510 5
[ bonds ]
C CT3
C +N
CT3 H1
CT3 H2
CT3 H3
O C
[ impropers ]
C CT3 +N O
[ CT3 ]
; this can also be done with the .c.tdb, but the atom naming is different
; and this can matter
[ atoms ]
N NH1 -0.470 0
HN H 0.310 1
CT CT3 -0.110 2
H1 HA 0.090 3
H2 HA 0.090 4
H3 HA 0.090 5
[ bonds ]
-C N
N HN
N CT
CT H1
CT H2
CT H3
[ impropers ]
N -C CT HN
; -C CT N -O ; commented to match with charmm psf
==============final rtp entries==================
==============final pdb======================
REMARK NONE *
REMARK DATE: 8/10/11 5:30:23 CREATED BY USER: ss2029
ATOM 1 CT3 ACE 1 -2.160 0.537 0.930 1.00 0.00 DIAL
ATOM 2 H1 ACE 1 -2.466 0.003 0.005 1.00 0.00 DIAL
ATOM 3 H2 ACE 1 -2.562 1.572 0.910 1.00 0.00 DIAL
ATOM 4 H3 ACE 1 -2.562 0.001 1.816 1.00 0.00 DIAL
ATOM 5 C ACE 1 -0.672 0.582 1.009 1.00 0.00 DIAL
ATOM 6 O ACE 1 -0.105 1.128 1.954 1.00 0.00 DIAL
ATOM 7 N ALA 2 0.000 0.000 0.000 1.00 0.00 DIAL
ATOM 8 HN ALA 2 -0.453 -0.444 -0.769 1.00 0.00 DIAL
ATOM 9 CA ALA 2 1.459 0.000 0.000 1.00 0.00 DIAL
ATOM 10 HA ALA 2 1.812 -0.497 0.897 1.00 0.00 DIAL
ATOM 11 CB ALA 2 1.949 -0.834 -1.207 1.00 0.00 DIAL
ATOM 12 HB1 ALA 2 1.514 -1.854 -1.155 1.00 0.00 DIAL
ATOM 13 HB2 ALA 2 1.628 -0.373 -2.167 1.00 0.00 DIAL
ATOM 14 HB3 ALA 2 3.056 -0.932 -1.214 1.00 0.00 DIAL
ATOM 15 C ALA 2 2.096 1.401 0.000 1.00 0.00 DIAL
ATOM 16 O ALA 2 1.425 2.432 0.000 1.00 0.00 DIAL
ATOM 17 N ALA 3 3.451 1.458 0.000 1.00 0.00 DIAL
ATOM 18 HN ALA 3 3.954 0.594 0.000 1.00 0.00 DIAL
ATOM 19 CA ALA 3 4.276 2.664 0.000 1.00 0.00 DIAL
ATOM 20 HA ALA 3 4.065 3.236 -0.897 1.00 0.00 DIAL
ATOM 21 CB ALA 3 3.864 3.539 1.207 1.00 0.00 DIAL
ATOM 22 HB1 ALA 3 2.776 3.756 1.155 1.00 0.00 DIAL
ATOM 23 HB2 ALA 3 4.063 3.014 2.167 1.00 0.00 DIAL
ATOM 24 HB3 ALA 3 4.407 4.508 1.214 1.00 0.00 DIAL
ATOM 25 C ALA 3 5.791 2.399 0.000 1.00 0.00 DIAL
ATOM 26 O ALA 3 6.597 3.328 0.000 1.00 0.00 DIAL
ATOM 27 N CT3 4 6.175 1.110 0.000 1.00 0.00 DIAL
ATOM 28 HN CT3 4 5.528 0.351 0.000 1.00 0.00 DIAL
ATOM 29 CT CT3 4 7.566 0.777 0.000 1.00 0.00 DIAL
ATOM 30 H1 CT3 4 8.048 1.201 0.907 1.00 0.00 DIAL
ATOM 31 H2 CT3 4 7.683 -0.327 0.000 1.00 0.00 DIAL
ATOM 32 H3 CT3 4 8.048 1.201 -0.907 1.00 0.00 DIAL
TER 33 CT3 4
END
==============final pdb======================
On Fri, Sep 9, 2011 at 11:35 AM, Mark Abraham <mark.abraham at anu.edu.au>wrote:
>
>
> On 09/09/11, *Sandeep Somani *<ssomani at gmail.com> wrote:
>
> Hi Mark
>
> Thanks for the rtp definitions. I added the rtp defs to aminoacids.rtp and
> renamed residue and atom names in pdb file.
>
> But pdb2gmx is still complaining about the terminal groups. It seems to
> have picked up ACE, but not CT3.
> pdb2gmx error message and modified pdb file are listed below. Any idea ? Do
> I need to modify .c.tdb ? If so, do you have that entry as well :) ?
>
>
> Ah, no. To get rid of the warning, you need to teach GROMACS that CT3 is a
> protein residue. (It knew about ACE from other force fields.)
> Copy residuetypes.dat from the share/top directory to your working
> directory, and add CT3 as Protein at the end.
>
> To get rid of the error, you need to tell pdb2gmx not to try to change your
> termini. Use pdb2gmx -ter, and choose "None" for both N and C termini.
>
> Mark
>
>
>
> Thanks
> Sandeep
>
> =============pdb2gmx error=======================================
> Back Off! I just backed up topol.top to ./#topol.top.3#
> Processing chain 1 (32 atoms, 4 residues)
> There are 3 donors and 3 acceptors
> There are 4 hydrogen bonds
> Identified residue ACE1 as a starting terminus.
> Warning: Residue CT32 in chain has different type (Other) from starting
> residue ACE1 (Protein).
> Identified residue ALA2 as a ending terminus.
> 8 out of 8 lines of specbond.dat converted successfully
> Start terminus ACE-1: NH3+
> End terminus ALA-2: COO-
>
> -------------------------------------------------------
> Program pdb2gmx_d, VERSION 4.5.4
> Source code file: pdb2top.c, line: 1070
>
> Fatal error:
> atom N not found in buiding block 1ACE while combining tdb and rtp
> For more information and tips for troubleshooting, please check the GROMACS
> website at http://www.gromacs.org/Documentation/Errors
> -------------------------------------------------------
> ================pdb2gmx error====================================
>
> ================modified pdb=====================================
> REMARK NONE *
>
> REMARK DATE: 8/10/11 5:30:23 CREATED BY USER: ss2029
>
> *ATOM 1 CT3 ACE 1 -2.160 0.537 0.930 1.00 0.00
> DIAL*
> *ATOM 2 HT1 ACE 1 -2.466 0.003 0.005 1.00 0.00
> DIAL*
> *ATOM 3 HT2 ACE 1 -2.562 1.572 0.910 1.00 0.00
> DIAL*
> *ATOM 4 HT3 ACE 1 -2.562 0.001 1.816 1.00 0.00
> DIAL*
> *ATOM 5 C ACE 1 -0.672 0.582 1.009 1.00 0.00
> DIAL*
> *ATOM 6 O ACE 1 -0.105 1.128 1.954 1.00 0.00
> DIAL*
> ATOM 7 N ALA 1 0.000 0.000 0.000 1.00 0.00
> DIAL
> ATOM 8 HN ALA 1 -0.453 -0.444 -0.769 1.00 0.00
> DIAL
> ATOM 9 CA ALA 1 1.459 0.000 0.000 1.00 0.00
> DIAL
> ATOM 10 HA ALA 1 1.812 -0.497 0.897 1.00 0.00
> DIAL
> ATOM 11 CB ALA 1 1.949 -0.834 -1.207 1.00 0.00
> DIAL
> ATOM 12 HB1 ALA 1 1.514 -1.854 -1.155 1.00 0.00
> DIAL
> ATOM 13 HB2 ALA 1 1.628 -0.373 -2.167 1.00 0.00
> DIAL
> ATOM 14 HB3 ALA 1 3.056 -0.932 -1.214 1.00 0.00
> DIAL
> ATOM 15 C ALA 1 2.096 1.401 0.000 1.00 0.00
> DIAL
> ATOM 16 O ALA 1 1.425 2.432 0.000 1.00 0.00
> DIAL
> ATOM 17 N ALA 2 3.451 1.458 0.000 1.00 0.00
> DIAL
> ATOM 18 HN ALA 2 3.954 0.594 0.000 1.00 0.00
> DIAL
> ATOM 19 CA ALA 2 4.276 2.664 0.000 1.00 0.00
> DIAL
> ATOM 20 HA ALA 2 4.065 3.236 -0.897 1.00 0.00
> DIAL
> ATOM 21 CB ALA 2 3.864 3.539 1.207 1.00 0.00
> DIAL
> ATOM 22 HB1 ALA 2 2.776 3.756 1.155 1.00 0.00
> DIAL
> ATOM 23 HB2 ALA 2 4.063 3.014 2.167 1.00 0.00
> DIAL
> ATOM 24 HB3 ALA 2 4.407 4.508 1.214 1.00 0.00
> DIAL
> ATOM 25 C ALA 2 5.791 2.399 0.000 1.00 0.00
> DIAL
> ATOM 26 O ALA 2 6.597 3.328 0.000 1.00 0.00
> DIAL
> *ATOM 27 N CT3 2 6.175 1.110 0.000 1.00 0.00
> DIAL*
> *ATOM 28 HN CT3 2 5.528 0.351 0.000 1.00 0.00
> DIAL*
> *ATOM 29 CT CT3 2 7.566 0.777 0.000 1.00 0.00
> DIAL*
> *ATOM 30 HT1 CT3 2 8.048 1.201 0.907 1.00 0.00
> DIAL*
> *ATOM 31 HT2 CT3 2 7.683 -0.327 0.000 1.00 0.00
> DIAL*
> *ATOM 32 HT3 CT3 2 8.048 1.201 -0.907 1.00 0.00
> DIAL*
> TER 33 CT3 2
> END
> ================modified pdb=====================================
>
>
>
>
> On Fri, Sep 9, 2011 at 2:52 AM, Mark Abraham <Mark.Abraham at anu.edu.au>wrote:
>
>> On 9/09/2011 8:05 AM, Sandeep Somani wrote:
>>
>>> Hi
>>>
>>> I am trying to set up a simulation for dialanine using charmm27 ff but am
>>> getting errors in pdb2gmx due to missing definitions in rtp file.
>>>
>>> I created the pdb file (below) using Charmm for the sequence
>>> ACE-ALA-ALA-CT3.
>>>
>>> pdb2gmx recognizes ALA atoms but not those of the terminal groups ACE and
>>> CT3.
>>>
>>
>> Yes, CHARMM in GROMACS has lacked these for some time. I asked over a year
>> ago for them to be added, but that hasn't happened. You can add
>>
>> [ ACE ]
>> [ atoms ]
>> CT3 CT3 -0.270 0
>> HT1 HA 0.090 1
>> HT2 HA 0.090 2
>> HT3 HA 0.090 3
>> C C 0.510 4
>> O O -0.510 5
>> [ bonds ]
>> C CT3
>> C +N
>> CT3 HT31
>> CT3 HT32
>> CT3 HT33
>> O C
>> [ impropers ]
>> C CT3 +N O
>>
>> [ CT3 ]
>> ; this can also be done with the .c.tdb, but the atom naming is different
>> ; and this can matter
>> [ atoms ]
>> N NH1 -0.470 0
>> HN H 0.310 1
>> CT CT3 -0.110 2
>> HT1 HA 0.090 3
>> HT2 HA 0.090 4
>> HT3 HA 0.090 5
>> [ bonds ]
>> -C N
>> N HN
>> N CAT
>> CT HT1
>> CT HT2
>> CT HT3
>>
>> [ impropers ]
>> N -C CAT HN
>> -C CAT N -O
>>
>> to the end of aminoacids.rtp in a local copy of the charmm27.ff folder to
>> make this work.
>>
>>
>>
>>> Does someone have the appropriate definitions for these residues ?
>>> ( I have tried http://swissparam.ch/ to generate an itp file for the
>>> molecule, but upon using it in GMX it does not give the same energy as that
>>> from charmm. )
>>>
>>
>> That will likely be a separate issue.
>>
>> Mark
>> --
>> gmx-users mailing list gmx-users at gromacs.org
>> http://lists.gromacs.org/**mailman/listinfo/gmx-users<http://lists.gromacs.org/mailman/listinfo/gmx-users>
>> Please search the archive at http://www.gromacs.org/**
>> Support/Mailing_Lists/Search<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-request at gromacs.org.
>> Can't post? Read http://www.gromacs.org/**Support/Mailing_Lists<http://www.gromacs.org/Support/Mailing_Lists>
>>
>
>
>
> --
> gmx-users mailing list gmx-users at 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-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20110910/8241c038/attachment.html>
More information about the gromacs.org_gmx-users
mailing list