[gmx-users] CHARMM36 force field available for GROMACS
Justin Lemkul
jalemkul at vt.edu
Mon Mar 10 03:26:56 CET 2014
On 3/9/14, 10:12 PM, Jim wrote:
> Thanks, but I'll have to ask you something again. :)
>
>
> On Sun, Mar 9, 2014 at 5:00 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>
>>
>>
>> On 3/9/14, 7:23 PM, Jim wrote:
>>
>>> Thanks for the clarification.
>>>
>>> *Bonds involving H5' and H5'' in RNA residues should indeed be changed
>>> to
>>> the H5'1 and H5'2 nomenclature.*
>>>
>>>
>>> Even though H5' and H5'' are connected to the same atom, what would be the
>>> problem to keep these H atom names H5' and H5'' instead of changing them
>>> to
>>> H5'1 and H5'2?
>>>
>>>
>> In principle, no, the CHARMM-specific names should work fine, but only if
>> you add corresponding entries in the merged.arn file to translate between
>> CHARMM and Gromacs nomenclature. We made the naming change to interface
>> smoothly with the Gromacs hydrogen-generation (.hdb) mechanism. If you
>> have a coordinate file that complies with CHARMM naming, i.e. from
>> CHARMM-GUI, then simple sed statements take care of any issues.
>
>
>
> Although this makes sense, it also brings other questions. I don't really
> understand the significance of the *hdb file - isn't its sole purpose to
> add hydrogens to coordinates that don't have all hydrogens? Second, what
Yes.
> kind of a file is merged.am? I've never seen such a file... Third, if the
It translates atom names between different nomenclature systems.
> H2' and H2'' of the RNA nucleotides do not correspond to the Gromacs
> nomenclature as well, then where exactly is their CHARMM->Gromacs
> conversion taken care of? In other words, if I use H5' and H5'' names
> instead of H5'1 and H5'2, isn't this equivalent to using the H2' and H2''
> names instead of H2'1 and H2'2 in terms of translating between
> nomenclatures?
>
Not technically. The difference is in connectivity. If one supplies pdb2gmx
with a coordinate file that has no hydrogens (i.e. a "typical" crystal
structure), atoms H5' and H5'' need to be generated on C5'. Gromacs handles
this using .hdb function type 6, and generates 2 hydrogens. There is no way to
assign names of H5' and H5'' unless one generates each separately, which in
theory, could be done but it's a bit awkward. If you specify atom H5' to be
generated as a base name, and two are generated, they are named H5'1 and H5'2.
In the case of H2' and H2'', they are not bonded to the same atoms. H2' is
bonded to C2' and thus can be generated as a single hydrogen very easily. H2''
is bonded to O2', so there is no problem with numerical suffices.
> Basically, I have a CHARMM 36 solvated file with all necessary hydrogens
> generated with CHARMM, do some rotation (for PBC matching) and then
> something like this:
>
> cat $pdbfile | sed 's/TIP3/SOL /' | sed 's/ OH2 SOL/ OW SOL/' | sed 's/
> H1 SOL/ HW1 SOL/' | sed 's/ H2 SOL/ HW2 SOL/' | sed 's/ SOD SOD/ NA NA /g
> ' | sed 's/ CLA CLA/ CL CL /g'
>
> I don't change the names of any other atom or residue. I use H5' and H5''.
> I change ADE, GUA, CYT, URA to use the H5' and H5'' atom names in the
> merge.rtp as well. I don't touch the [bonds] because they already use these
> atom names. I generate a topology and *.gro with pdb2gmx. It looks like
> this (e.g. for one residue):
>
> [ atoms ]
>
> ; nr type resnr residue atom cgnr charge mass typeB
> chargeB massB
>
> ; residue 17 GUA rtp GUA q -0.5
>
> 1 ON5 17 GUA O5' 1 -0.66 15.9994 ;
> qtot -0.66
>
> 2 HN5 17 GUA H5T 2 0.43 1.008 ;
> qtot -0.23
>
> 3 CN8B 17 GUA C5' 3 0.05 12.011 ;
> qtot -0.18
>
> 4 HN8 17 GUA H5' 4 0.09 1.008 ;
> qtot -0.09
>
> 5 HN8 17 GUA H5'' 5 0.09 1.008 ;
> qtot 0
>
> 6 CN7 17 GUA C4' 6 0.16 12.011 ;
> qtot 0.16
>
> 7 HN7 17 GUA H4' 7 0.09 1.008 ;
> qtot 0.25
>
> 8 ON6B 17 GUA O4' 8 -0.5 15.9994 ;
> qtot -0.25
>
> 9 CN7B 17 GUA C1' 9 0.16 12.011 ;
> qtot -0.09
>
> 10 HN7 17 GUA H1' 10 0.09 1.008 ;
> qtot 0
>
> 11 NN2B 17 GUA N9 11 -0.02 14.007 ;
> qtot -0.02
>
> 12 CN5 17 GUA C4 12 0.26 12.011 ;
> qtot 0.24
>
> 13 NN1 17 GUA N2 13 -0.68 14.007 ;
> qtot -0.44
>
> 14 HN1 17 GUA H21 14 0.32 1.008 ;
> qtot -0.12
>
> 15 HN1 17 GUA H22 15 0.35 1.008 ;
> qtot 0.23
>
> 16 NN3G 17 GUA N3 16 -0.74 14.007 ;
> qtot -0.51
>
> 17 CN2 17 GUA C2 17 0.75 12.011 ;
> qtot 0.24
>
> 18 NN2G 17 GUA N1 18 -0.34 14.007 ;
> qtot -0.1
>
> 19 HN2 17 GUA H1 19 0.26 1.008 ;
> qtot 0.16
>
> 20 CN1 17 GUA C6 20 0.54 12.011 ;
> qtot 0.7
>
> 21 ON1 17 GUA O6 21 -0.51 15.9994 ;
> qtot 0.19
>
> 22 CN5G 17 GUA C5 22 0 12.011 ;
> qtot 0.19
>
> 23 NN4 17 GUA N7 23 -0.6 14.007 ;
> qtot -0.41
>
> 24 CN4 17 GUA C8 24 0.25 12.011 ;
> qtot -0.16
>
> 25 HN3 17 GUA H8 25 0.16 1.008 ;
> qtot 0
>
> 26 CN7B 17 GUA C2' 26 0.14 12.011 ;
> qtot 0.14
>
> 27 HN7 17 GUA H2'' 27 0.09 1.008 ;
> qtot 0.23
>
> 28 ON5 17 GUA O2' 28 -0.66 15.9994 ;
> qtot -0.43
>
> 29 HN5 17 GUA H2' 29 0.43 1.008 ;
> qtot 0
>
> 30 CN7 17 GUA C3' 30 0.01 12.011 ;
> qtot 0.01
>
> 31 HN7 17 GUA H3' 31 0.09 1.008 ;
> qtot 0.1
>
> 32 ON2 17 GUA O3' 32 -0.57 15.9994 ;
> qtot -0.47
>
> ; residue 18 GUA rtp GUA q -1.0
>
>
> I use that topology and the GRO file to run MD. Is there anything incorrect
> here or do I have to have done anything additional, i.e. deal with some
> merge.am naming conversion?
>
If no fatal errors came up, the result is fine. The .arn mechanism is just for
convenience when running pdb2gmx.
We will be releasing a new version of the force field in a day or two that
corrects the errors you have noted, expands the .hdb mechanism to be fully
integrated with RNA, and fixes a few small errors (nothing serious).
-Justin
--
==================================================
Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 601
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul
==================================================
More information about the gromacs.org_gmx-users
mailing list