Moltemplate OPLS-AA improper parameters

LAMMPS LAMMPS Mailing List Mirror
1

Firstly, I would like you to thank you for making the Moltemplate tool
available. It has been invaluable to me. The reason I'm writing is because
of what appears to me to be a couple of errors in the specification of
improper torsions for the OPLS-AA force field in oplsaa.lt, but perhaps I am
misunderstanding something.

Firstly, the "harmonic" improper_style is used, when the OPLS-AA force field
uses the same type of function for impropers as the AMBER force field, which
is equivalent to the "cvff" improper_style in LAMMPS.

Jason Lambert wrote the conversion script. ("tinkerparm2lt.py") I
don't know why Jason chose "harmonic" instead of "cvff" but I strongly
suspect it is because improper_style "harmonic" can be used with GPU
acceleration (harmonic/kk), and "cvff" can not (yet). Since, for small
angles the "harmonic" function is a pretty good approximation of the
cosine function used in "improper_style cvff", this is a pretty
important reason to keep things as they are.
    If the entire purpose of the improper interaction is to keep the 4
atoms coplanar (which it is), then hopefully it does not matter too
much if you use one or the other (as long as I did not screw up the
magnitude of the spring constants. I think they are correct, but it's
possible they might be off by a factor of 2. Let me know if you spot
something.) I should also admit that we also ignore the
"Urey-Bradley" terms in the OPLS force field. (Although I think there
is only one of them, and apparently most people don't bother with
these terms unless they care about the high frequency vibrational
spectra.)

   Nevertheless, if you have reason to believe that it's important to
use a cosine function instead of a quadratic one, I can make an
alternative version of the oplsaa.lt file which uses "cvff" instead.
Nobody has asked for this so far, but perhaps they are just not as
diligent as you. I'm happy to do this.

Secondly, the
parameters for the improper torsions do not seem to match those in the
papers in which the OPLS-AA force field is specified - in particular, the
improper torsion parameters are apparently from the original AMBER paper,
and the latest revision of OPLS (OPLS3) (published in
http://dx.doi.org/10.1021/acs.jctc.5b00864) appears to use similar
parameters. I realise you have taken these parameters from the TINKER
parameter file for OPLS, but they do not seem to be correct.

The OPLSAA parameters in the "oplsaa.lt" file were taken from this file:
https://dasher.wustl.edu/tinker/distribution/params/oplsaa.prm
Here's an excerpt from the most recent version of that file:
"The parameters supplied with TINKER are from "OPLS All-Atom Parameters
for Organic Molecules, Ions, Peptides & Nucleic Acids, July 2008" as
provided by W. L. Jorgensen, Yale University during June 2009. These
parameters are taken from those distributed with BOSS Version 4.8."

It's really confusing. There are multiple sets of OPLSAA parameters,
as they have changed over the years. I am under the impression that
they update the parameters frequently and they don't always publish a
new paper when they do this (but perhaps I'm wrong). The parameters
have definitely changed a lot from the original OPLS papers.

I don't know if this helps, but here is a link to the parameters used
in BOSS 4.9. Hopefully they are similar to BOSS 4.8. The list of
parameters is located on pages 66-74 of this file:
http://zarbi.chem.yale.edu/doc/boss49.pdf

Unfortunately, I could not tell if there are any _improper_ torsion
parameters in this file (or just ordinary torsion parameters).

If you can find a recent source for OPLSAA parameters which disagrees
with the moltemplate parameters in "oplsaa.lt" file, let me know and I
will definitely fix it. There are only 5 improper interactions in
this version of OPLSAA. If that's all that's wrong, it is easy to fix
them manually.

Thanks for the kindly words about moltemplate.

Cheers

Andrew

2

The OPLS-AA parameters in Gromacs are provided in the share/top/oplsaa.ff
directory of the gromacs code distribution, which can be downloaded from
http://www.gromacs.org/Downloads
The bonded parameters are in ffbonded.itp and the nonbonded parameters
in ffnonbonded.itp. Gromacs uses units of kJ/mol for energies and nm for
distances, so there are differences in the parameters to those in oplsaa.lt
due to units. Also, the harmonic bond length and angle force constant
parameters differ by a factor of two from those for the LAMMPS "harmonic
style" and the proper dihedrals use the Ryckaert-Bellemans function, so
the dihedral parameters in Gromacs are related to linear combinations of
those for the LAMMPS "opls" dihedral style (the Gromacs manual gives
the formulae relating the OPLS and Ryckaert-Bellemans parameters).

For the systems I am studying, I only spotted inconsistencies between
the OPLS papers and Gromacs parameters vs moltemplate parameters for the improper dihedrals.

Grumbling because I've never been able to locate a complete list of
force field parameters in any of the OPLS papers or their supplemental
files. (It seems to be common practice. Perhaps they just like to
announce that new software tools are available. I do that too.)
Anyway, thanks for the reference to the GROMACS
"oplsaa.ff/ffbonded.itp" file. That helps. Here's the link:
https://github.com/gromacs/gromacs/blob/master/share/top/oplsaa.ff/ffbonded.itp
For reference, here's the TINKER "oplsaa.prm" file we get our parameters:
https://dasher.wustl.edu/tinker/distribution/params/oplsaa.prm

Alas, now I am more confused than before, since there are a number of
disagreements between these files. I'll make some guesses about where
the these files agree and the two files differ.

1) AGREEMENT: This line from "oplsaa.ff/ffbonded.itp":
define improper_O_C_X_Y 180.0 43.93200 #(gromacs)
     ...seems to correspond to these lines from moltemplate and tinker:
improper_coeff @improper:X-X-3-4 harmonic 10.5 180.0 #(moltemplate)
imptors 0 0 3 4 21.000 180.0 2 #(tinker)
     It also seems to corresponds to these lines as well:
improper_coeff @improper:X-X-3-52 harmonic 10.5 180.0 #(moltemplate)
imptors 0 0 3 52 21.000 180.0 2 #(tinker)

     (Notes to self: A) gromacs files use kJ units instead of kCal.
B) Improper energies in tinker files are larger because they do not
absorb the factor of 2 into the coefficient. C) The central "hub"
atom is the 2nd atom in the list in gromacs files, but it is the 3rd
atom in the list in tinker and the corresponding moltemplate files.)

2) CONFLICT: This line from "oplsaa.ff/ffbonded.itp":
define improper_Z_N_X_Y 180.0 4.18400 #(gromacs)
   ...seems to disagree with these lines from moltemplate and tinker:
improper_coeff @improper:X-X-24-X harmonic 2.5 180.0 #(moltemplate)
imptors 0 0 24 0 5.000 180.0 2 #(tinker)
   If the gromacs file is correct, I should change the 2.5 to 1.0

3) AGREEMENT: This line from "oplsaa.ff/ffbonded.itp":
define improper_Z_CM_X_Y 180.0 62.76000 2 #(gromacs)
     ...seems to correspond to these lines from moltemplate and tinker:
improper_coeff @improper:X-X-47-X harmonic 15.0 180.0 #(moltemplate)
imptors 0 0 47 0 30.000 180.0 2" #(tinker)

4) CONFLICT: This line from "oplsaa.ff/ffbonded.itp":
define improper_Z_CA_X_Y 180.0 4.60240 2 #(gromacs)
   ...seems to disagree with these lines from moltemplate and tinker:
improper_coeff @improper:X-X-48-X harmonic 2.5 180.0 #(moltemplate)
imptors 0 0 48 0 5.000 180.0 2 #(tinker)
   If the gromacs file is correct, I should change the 2.5 to 1.1

5) CONFLICT: The next line from "oplsaa.ff/ffbonded.itp":
define improper_N2_X_N2_N2 180.0 43.93200 2 #(gromacs)
    ...does not correspond to anything in the gromacs or tinker files.
    David, should I create a new line in the moltemplate file like this?:
improper_coeff @improper:24_24_X_24 harmonic 10.5 180.0

    (Background information: I chose "24" because it corresponds to
"i24" which appears to correspond to several different types of
nitrogen atoms (all of which are hopefully equivalent to "N2" in the
gromacs file) according to these equivalence rules (which originally
come from the "oplsaa.prm" and replace atom types 179,180,181... with
"i24":
set type @atom:179 charge -0.76 # "Amide -CO-NH2"
set type @atom:180 charge -0.5 # "Amide -CO-NHR"
set type @atom:181 charge -0.14 # "Amide -CO-NR2"
replace{ @atom:179 @atom:179_b24_a24_d24_i24 }
replace{ @atom:180 @atom:180_b24_a24_d24_i24 }
replace{ @atom:181 @atom:181_b24_a24_d24_i24 }
    Let me know if you think I guessed the atom types incorrectly or
if I am misunderstanding this entirely. (My chemistry may be rusty,
but it seems weird to have 3 nitrogen atoms connected to the same
atom. Perhaps this does not happen often. Should we just ignore this
weird rule?)

6) CONFLICT: The next line from "oplsaa.ff/ffbonded.itp":
define improper_X_NO_ON_NO 180.0 43.93200
    ...does not seem to correspond to anything in the gromacs or tinker files.
   I am not certain what "NO" and "ON" correspond to. Again, in
gromacs files, the second atom ("NO") is the central hub in the
interaction. In TINKER(oplsaa.prm), and hence, in
MOLTEMPLATE(oplsaa.lt), it looks like these might be the relevant atom
types:
set type @atom:701 charge 0.54 # "Nitroalkane -NO2" (NITROGEN)
set type @atom:702 charge -0.37 # "Nitroalkane -NO2" (OXYGEN)
replace{ @atom:701 @atom:701_b102_a102_d102_i102 }
replace{ @atom:702 @atom:702_b103_a103_d103_i103 }
    Hence I perhaps I could add a rule like this to moltemplate:
improper_coeff @improper:X_103_102_102 harmonic 10.5 180.0
    ...but I won't add this rule unless someone can reassure me that I
guessed the atom types correctly (Note to self: remember we have to
swap the 2nd and 3rd atoms in the list, because gromacs and tinker
disagree on atom order)

   I'm happy to make changes #2 and #4. I'm curious to hear your
opinion #5 and #6.

Cheers, and thanks again for reporting the differences.

Andrew