This sounds like a bug in my code from long ago, but I don't remember
enough to comment in detail. One thing to make sure of is that
impropers aren't being overlooked because of bad comparisons between
molecule and FF atom types. It may be that it's just an issue with
ordering the atoms in the prmtop, which would be less serious than
missing terms.
Bill
Reinis Danne <gam4rei.gmail.com> wrote:
> Hi!
>
> I'm building saccharide monomer 0BC from Glycam 06h and found
> that tleap does not put atoms of improper dihedral in correct
> order. Moreover it isn't even consistent within one molecule
> where atoms of the same improper type get ordered differently.
>
> According to AmberTools and PREP file documentation atoms of
> impropers are ordered as follows:
> - Central atom is always in 3rd position.
> - The three atoms bonded to the central atom are ordered by
> their atom type.
> - In case of identical atom types, they are sorted by atom
> index.
>
> But tleap produces this:
> B 65: 1.500 3.14 2.0 :2.C2N :2.H2N :2.N2 :2.C2 (30,29,28,5)
> B 66: 1.500 3.14 2.0 :2.C4N :2.C4 :2.N4 :2.H4N (20,9,18,19)
>
> With Glycam 06g it put both of them in the same order:
> B 65: 1.500 3.14 2.0 :2.C2 :2.C2N :2.N2 :2.H2N (26,30,28,29)
> B 66: 1.500 3.14 2.0 :2.C4 :2.C4N :2.N4 :2.H4N (12,16,14,15)
>
> The difference is that Glycam 06g had this improper specified
> using wildcard atoms (X X N H) while in version 06h it is given
> explicitly (C H Ng Cg).
>
> So tleap actually has another rule for atom ordering in impropers:
> - Atoms matched by wildcards are moved at the beginning of the
> list and sorted by index.
>
> Is it intentional and as intended? If so documentation should be
> amended.
>
> I think that atoms should be ordered by their type and index as
> tie-breaker regardless if it happened to be matched by improper
> with wildcards in it or not. That would be consistent with
> documentation, but tleap should be corrected.
>
> In any case, there is an inconsistency somewhere which results
> in different order of the same improper type in one molecule.
>
> The correct order would be this:
> B 65: 1.500 3.14 2.0 :2.C2N :2.C2 :2.N2 :2.H2N (30,5,28,29)
> B 66: 1.500 3.14 2.0 :2.C4N :2.C4 :2.N4 :2.H4N (20,9,18,19)
>
> Or if ordered as in parameter file:
> B 65: 1.500 3.14 2.0 :2.C2N :2.H2N :2.N2 :2.C2 (30,29,28,5)
> B 66: 1.500 3.14 2.0 :2.C4N :2.H4N :2.N4 :2.C4 (20,19,18,9)
>
> The cause of different order within the same molecule in case of
> Glycam 06h probably is a problem in
> leap/parmSet.c:ParmSetImproperOrderAtoms() or in the loop in
> leap/unitio.c:zbUnitIOIndexTorsionParameters() from which it is
> called.
>
> Here is the explanation as I understand it:
> - Improper parameters in leap are stored in TORSIONPARMt struct.
>
> - This struct has 'sOrder' field.
>
> - sOrder is basically a map between atom types and atom indexes
> in the improper or in case of the parameters it maps sorted
> atom types to the order they were in the original parameter
> file.
>
> - When improper parameters are read, their atom types are sorted
> according to the rules above and sOrder is changed to reflect
> how the order changed.
>
> - When improper parameters are added to specific torsion,
> atom types of the relevant atoms are sorted according to the
> rules above and sOrder is changed. Since only the order of
> atom types was changed, but not the order of the atoms, sOrder
> now tells which atom types correspond to which atoms. sOrder
> is stored in the torsion with the new atom type order and old
> atom order.
>
> - After that ParmSetImproperOrderAtoms() is called and should
> change atom order to correspond to the order of atom types as
> shown by sOrder.
>
> - In this function the atoms which are given explicitly (without
> wildcards) are ordered first without looking at sOrder at all.
>
> - Then both, atoms and their types, are reordered according to
> sOrder.
>
> - After that atoms matching wildcards are moved to the beginning
> and sorted by their index.
>
> - And finally atoms with identical atom types are sorted by
> their index.
>
> - After this in unitio.c:zbUnitIOIndexTorsionParameters() the
> function iParmSetAddImproperTerm() is called which resets
> sOrder to '0123' since now it receives already sorted atom
> types. In the next loop atom order will not be changed to the
> one in original parameter file.
>
> So the issue is that sOrder is changed in between the two calls
> to ParmSetImproperOrderAtoms() for the same improper. For some
> reason iParmSetAddImproperTerm() is called again (the first time
> was when loading parameter file) and now it finds sOrder='0123'
> when in the first case it was '0321'. So atoms of the second
> improper are not reordered leading to different atom orders for
> the two impropers.
>
> I'm not familiar with inner workings of leap so I don't know
> what would be the correct fix for this.
>
> I have not checked if there is similar issue to proper dihedrals
> and other terms.
>
> Any comments and suggestions would be appreciated. Also,
> pointing out if I got something wrong would be helpful.
>
>
> Best Regards,
> Reinis Danne
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed May 08 2013 - 18:00:02 PDT
