I have been having a hard time figuring out where to get the OPLS dihedral
parameters for lammps. I see in the OPLS dihedral form used by lammps, one
specifies the V1 V2 V3 and V4 in the dihedral_coeff command. However, in all
the OPLS parameter sets I have looked at, the one given in tinker being what
I think is close to the latest version of OPLS
(http://dasher.wustl.edu/ffe/distribution/params/oplsaa.prm), I only see a
mention of what I believe to be V1 V2 V3 and f1 f2 f3 where the latter are
the phase angle offsets.
...
or is there some
conversion between the V1 V2 V3 f1 f2 f3 to the V1 V2 V3 V4 sets that I am
missing.
You are asking for details about the format of the "oplsaa.prm" file,
which is is not part of LAMMPS. It is a file distributed with the
TINKER software, and used by MOLTEMPLATE.
Since I wrote moltemplate, I will do my best to answer.
If you discover that my answer is incorrect, please let me know.
However, this file is distributed by the Ponder lab, and they are the
authoritative source for information about this file.
Anyway, in the 2001 JPCA paper cited by the LAMMPS manual,
(http://courses.chem.psu.edu/chem408/reading/MM_topics/jorgensen_perfluoro_jpca_2001.pdf)
the formula for the torsion energy is:
0.5*(V1*(1+cos(x+f1) + V2*(1-cos(2x+f2)) + V3*(1+cos(3x+f3)) +
V4*(1-cos(4x+f4)))
There is no general way to convert this formula to the formula used by
LAMMPS' "dihedral_style opls", for arbitrary values of f1,f2,f3,f4
http://lammps.sandia.gov/doc/dihedral_opls.html
For dihedral style opls, the energy is:
0.5*(K1*(1+cos(x)) + K2*(1+cos(2x)) + K3*(1+cos(3x)) + K4*(1+cos(4x)))
(Note there is a sign change in the two formulas.: "1-cos(2x+f2)"
instead of "1+cos(2x)".)
However, it turns out that, for every torsion interaction in OPLSAA:
f1=0
f2=180 # (this negates the sign change)
f3=0
f4=180
(See the torsion section of that "oplsaa.prm" file.) Consequently,
the "dihedral_style_opls" formula is general enough, and the
conversion you are looking for is:
K1 = V1
K2 = V2
K3 = V3
K4 = V4
(This is what Jason Lambert's "oplsaa_moltemplate.py" script does when
it converts the "oplsaa.prm" file into moltemplate format
(oplsaa.lt).)
Finally, I am having a hard time interpreting the parameters in the torsional part of
the oplsaa file. For instance,
torsion 13 3 20 13 4.669 0.0 1 5.124 180.0 2 0.000 0.0 3
I am assuming 4.669,5.124,and 0.000 correspond to a V1 V2 V3. Is this true?
Yes.
What are the 1 2 3 referring to along with the 0.0 180.0 0.0 (I am assuming these are a f1,f2,f3).
Yes.
(The "1", "2" and "3", are the integer frequencies in the Fourier
expansion. They never vary in the "oplsaa.prm" file. They are always
1,2,3.)
Unfortunately, the "oplsaa.prm" file does not include V4 parameters.
However, they are not needed if you are using the original, basic
OPLSAA force-field. In that case V4=0. (Apparently V4 is non-zero
for the perfluouroalkanes discussed in the 2001 paper, but not in the
original 1996 OPLSAA paper. See below.)
The Watkins paper cited in the lammps manual only
provides parameters for perfluouroalkanes, so I was wondering where I need
to go to obtain these parameters for other molecules,
The original OPLS paper is from 1996, and it is not only for perfluouroalkanes
http://pubs.acs.org/doi/abs/10.1021/ja9621760
It uses the same formula for the dihedral energy.
Incidentally, you do not have to use "dihedral_style opls" to simulate
molecules using the OPLS force-fields in LAMMPS. If you want a more
general formula, I recommend "dihedral_style fourier".
http://lammps.sandia.gov/doc/dihedral_fourier.html
That dihedral style will allow you to customize the f1,f2,f3,f4,...
offsets, and the frequencies as well (1,2,3,4,...)
If you need more detail than this, the TINKER mailing list is probably
the best source of information.
Andrew