[gmx-users] Cannot Generate CHARMM36 Topology for POPC with Virtual Sites

[Justin Lemkul]

On 4/8/14, 12:45 PM, Chris Ing wrote:
> I'm trying to build a charmm36.ff topology file for a POPC with virtual
> sites. In GROMACS 4.5.5 and 4.6.3 (with the most recent CHARMM36 March
> 2014), when I run:
>
> pdb2gmx -f POPC.pdb -o POPC.gro -p POPC.top -i POPC.posre -vsite hydrogen
> grompp -f empty.mdp -c POPC.gro -p POPC.top -o test.tpr
>
> I get multiple constraint errors:
>
> ERROR 1 [file POPC.top, line 311]:
>
>    No default Constr. No Conn. types
>
> Which is due to the lack of two constraints not in ffbonded.itp:
>
> MCH3 NTL
> MCH3 CTL2
>
> MCH3 is the virtual site of a CH3 group, and NTL refers to the choline
> nitrogen in the POPC headgroup and CTL2, in this case, is the carbon
> closest to the CH3 group on the POPC tail.
>
> The same question was asked in the gmx-users thread "amber99sb_virtual site
> lipids dopc MCH3_N constraint types" by Song Ke, where the reply from Erik
> Marklund was to manually add this value "from bond lengths, center of mass
> for the group and the moment of inertia.". In lieu of wiping an inch of
> dust off my classical mechanics textbook, I wondered if I could copy some
> additional constraint values from the OPLS ffbonded.itp file.
>
> The charmm36.ff ffbonded.itp file states that the constrainttypes section
> is copied directly from OPLS (the restraint values are slightly adjusted?).
> Yet, many entries were not copied from OPLS to CHARMM36:
>
> OPLS:
>
> MNH3   CT_2 2    0.158255
> MNH3   MNH3 2    0.080236
> MCH3A     C 2    0.166040
> MCH3A   C_2 2    0.166040
> MCH3A    CW 2    0.164312
> MCH3A    CV 2    0.164312
> MCH3A    CS 2    0.163448
> MCH3A    CA 2    0.164888
> MCH3A    CB 2    0.164888
> MCH3A     N 2    0.159040
> MCH3A     S 2    0.193874
> MCH3A MCH3A 2    0.092160
> MCH3B    CT 2    0.167031
> MCH3B  CT_2 2    0.167031
> MCH3B  CT_3 2    0.167031
> MCH3B    CO 2    0.167031
> MCH3B MCH3B 2    0.091456
>
> CHARMM36:
>
> MNH3     CT3   2    0.159603
> MNH3    MNH3   2    0.082619
> MCH3     CT1   2    0.168122
> MCH3     CT2   2    0.167162
> MCH3     CT3   2    0.167354
> MCH3       S   2    0.195314
> MCH3    MCH3   2    0.093582
>
> Would it be acceptable to use these OPLS constraint values in CHARMM36 to
> resolve my problem?:
>
> MCH3A     N 2    0.159040
> MCH3B    CT 2    0.167031
>
> Thanks gmx-users/Justin Lemkul!
>

For what it's worth, we have done exactly zero testing with virtual sites.  The 
sections of text you are referring to were copied basically verbatim from what 
was distributed as charmm27.ff in Gromacs prior to our creation of charmm36.ff. 
  I don't know if that's of any use, but please consider it at least a strong 
disclaimer :)

> -Chris
>
> P.S. I had similar missing dihedral/ub errors when trying to build my
> charmm36.ff protein topology with virtual sites whenever I used a neutral
> terminus. I didn't look into the details because using a charged terminus,
> I had no errors.
>

If you can provide me with specific errors and the atomtypes to which they 
correspond, I will look into it.  We recently fixed a problem with missing U-B 
types for neutral termini with branched-chain amino acids, but the necessary 
parameters were added before our March release.  To my knowledge, that issue 
should be resolved, but if there are remaining problems, I'd be happy to look 
into it.  If it only occurs when using virtual sites, again, that's because we 
never looked into it.

-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

==================================================