On Wed, Jun 5, 2013 at 3:05 AM, FyD <fyd.q4md-forcefieldtools.org> wrote:
> Dear All,
>
> If one loads the same number of _improper_ dihedrals (FF parameters)
> in LEaP using either the generic (i.e. improper defined with 'X'
> letters) or the specific 'way' (no 'X' letter) LEaP does not
> rigorously apply the same impropers...
>
> I take an example:
> - I first executed leaprc.ff99SB (so there are generic improper
> dihedrals in parm99.dat) and saved the prmtop/prmcrd files for
> ACE-TYR-NME: 10 impropers are generated.
>
> - Then I loaded my own frcmod file with only specific impropers; see below:
> CA-CA-C -OH 1.10000000e+00 180.0 2. taken from parm99.dat
> CA-CA-CA-CT 1.10000000e+00 180.0 2. taken from parm99.dat
> C -CT-N -H 1.10000000e+00 180.0 2. taken from parm99.dat
> CT-N -C -O 1.05000000e+01 180.0 2. adapted from
> parm99.dat
> CA-CA-CA-HA 1.10000000e+00 180.0 2. adapted from
> parm99.dat
> C -CA-CA-HA 1.10000000e+00 180.0 2. adapted from
> parm99.dat
> and saved the prmtop/prmcrd files for the same ACE-TYR-NME dipeptides
>
> I compare the generated prmtop files:
>
> diff generic.top specific.top
> 230c230
> < 39 48 -42 -45 18 54 42 -48 -51 18
> ---
> > 39 48 -42 -45 18 42 54 -48 -51 18
> -> One observes differences in the %FLAG DIHEDRALS_INC_HYDROGEN
> section of the prmtop file
>
> Then, I used 'rdparm' to study these 10 impropers; saved the 'rdparm'
> output in two text files, that I compared:
>
> diff generic.txt specific.txt
> 63c63
> < B 62: 1.100 3.14 2.0 :1.CZ :1.CD1 :1.CE1 :1.HE1 (19,15,17,18)
> ---
> > B 62: 1.100 3.14 2.0 :1.CD1 :1.CZ :1.CE1 :1.HE1 (15,19,17,18)
>
> A single CA-CA-CA-HA improper is different; however, there are 4
> CA-CA-CA-HA impropers in the Tyrosine residue; this means 3 are
> applied using an homogeneous way, and for one the two first atoms are
> inverted.
>
> Obviously I understand that the impact on the energy value is
> negligible, but this is quite painful when one tries to generate
> identical prmtop files:
>
> Any idea where/how to correct that in the LEaP source code?
>
I had noticed that myself not too long ago. In 'real' impropers, the only
atom whose position 'matters' is the third (central) atom -- the order of
the others is arbitrary. However, since we treat impropers the same way as
we treat proper dihedrals, the overall order will impact the energy a
little bit. Therefore, the convention is applied to order the atoms in
alphabetical order by type.
However, there is no rigorously correct way of sorting 2 identical types
(CA in your case). As a result, there is no well-defined ordering in this
particular case (note that one of the CAs will always be the third atom, it
is only the other 2 that can swap).
A more rigorous way of fixing the ordering is to order based on alphabetic
order of types (to maintain 'classical' behavior of all force fields), then
order degenerate types based on their atom names. While it is fairly
arbitrary, it will at least yield reliably consistent topology files. And
as you mentioned, improper ordering will have a minor effect on energies,
but an even smaller effect on dynamics, I'd imagine.
A good place to start may be the zParmSetOrderImproperAtoms function
(parmSet.c).
HTH,
Jason
--
Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
352-392-4032
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Jun 05 2013 - 07:30:03 PDT
