[gmx-users] how likely that charmm-in-gromacs LJ-14 interactions are affected by inclusion of cgenff stuff when adding a new molecule with published charmm parameters?

[Justin Lemkul]

On 8/8/16 1:05 AM, Christopher Neale wrote:
>
> Dear Justin:
>
> I have another question, this time about LJ 1-4 interactions. I notice that
> charmm36-jun2015.ff/ffbonded.itp has the following in [ pairtypes ]
>
> NH1     S  1  0.316269044940  1.255200000000
>
> and it has the following in [ atomtypes ]
>
> S    16    32.060000    0.000  A  0.356359487256  1.88280 NH1     7
> 14.007000    0.000  A  0.329632525712  0.83680
>
> charmm in gromacs is comb-rule 2, so sigma = aritmetic and epsilon =
> geometric
>
> For the S-NH1 LJ1-4 interaction, the epsilon makes sense, but the sigma is
> off.
>
> Am I to conclude from this that charmm uses a whole bunch of semi-arbitrary
> LJ-14 terms? I checked in gromacs by simply removing the [ pairtypes ]
> section and when I did that my energy match from gromacs to charmm executable
> goes from not bad to really terrible (which aligns with the idea of charmm
> having a bunch of pre-set 1-4 interactions).
>
> My confusion stems from m reading of charmm mailing list posts in which alex
> states that modifications of the LJ1-4 scaling is only done in special cases
> (e.g,
> https://www.charmm.org/ubbthreads/ubbthreads.php?ubb=showflat&Number=1048 )
>

There are a few specific interactions that use special LJ 1-4 parameters.  What
is listed above is correct, because NH1 has a different Rmin/2 when in 1-4
interactions.  These parameters are explicitly listed in par_all36_prot.prm as
extra columns, as explained in the CHARMM documentation.  So these interactions
aren't "arbitrary" or hidden from the user, but you do have to know what you're
looking at.

> I ask because I am trying to ensure that the inclusion of all the cgenff
> stuff in the charmm36-jun2015.ff gromacs distribution is not going to
> unexpectedly affect LJ1-4 interactions in a molecular modification for which
> parameters have been published by another group who does their simulations
> with charmm. In this case, my charmm vs gromacs executable LJ difference goes
> up from ~ 0.01 kcal/mol for a 200 aa protein to 0.08 kcal/mol simply by
> adding a ~10 heavy atom covalent acylation.
>

CGenFF and the parent CHARMM force field should be independent.  One should 
*not* be using CGenFF atom types for constituent residues of a protein, only 
noncovalent interactions like ligands, co-solvents, etc.  So CGenFF certainly 
will not affect 1-4 interactions in the protein.

You can verify the 1-4 vs. SR interactions in CHARMM with:

! list all 1-4 interactions and energies
! this actually returns 2*energy
coor dist ener cut 999. nononbonds noexclusions 14exclusions sele all end sele 
all end

! just the short-range energies
! likewise, 2*energy
coor dist ener cut 999. nonbonds noexclusions no14exclusions sele all end sele 
all end

This sometimes points to the nature of the underlying problem.

Without knowing exactly what you're trying to reproduce (what paper, what 
interactions, etc) there's not much more I can say.  Beware - some people play 
weird tricks that "work" but are not necessarily sensible.  This is why we keep 
very tight control over what actually gets integrated into the CHARMM force 
field.  The SI of some papers scares us...

-Justin

-- 
==================================================

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

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