AMBER force field dihedral in LAMMPS

Dear all,

When I read the LAMMPS manual about the dihedral style CHARMM, I came up with a question. The manual says the 1-4 term in CHARMM style dihedral has to be zero when using AMBER force field. If I understand correctly, when this scale factor is set to zero, 1-4 non-bonded terms in dihedrals are firstly turned off and then evaluated by the "special_bonds amber” command. The "special_bonds amber” globally evaluates all the 1-4 non-bonded term for AMBER dihedral. However, when I went through the definition of topology file in AMBER, the instruction says that, for some particular dihedrals, their 1-4 non bonded terms are turn off. Basically, it means some of the 1-4 terms in dihedrals are evaluated by 0.5 and 0.8333, while some are skipped. I am wondering how can I achieve this in LAMMPS.

Any suggestion will be appreciated.
Andy

no. you are not understanding this correctly.

the special_bonds command does not do any computation, it merely sets the scaling factor for 1-2, 1-3, and 1-4 non-bonded interactions as computed by the pair style and controls how atoms are detected at 1-2, 1-3, or 1-4 pairs.

the charmm dihedral style is special for historical reasons. it predates a more capable implementation of the detection of special_bonds atom pairs and has not been updated since that is a quite complex process due to the complexity of the CHARMM force field.
thus to work around the limitation of the older LAMMPS versions, the 1-4 pair style interactions were excluded from the pair style and computed as part of the dihedral style instead. this requires the scaling factor of 1.0, except when to handle double counting, e.g. in opposing benzene ring atoms, when it needs to be set to 0.5.

since with amber force fields, you do regular processing of 1-4 interactions as part of the pair style, you must set the factor in the charmm dihedral style to 0.0 so it doesn’t try to compute the 1-4 interactions a second time, or rather fail because the CHARMM 1-4 treatment in the dihedral style also requires a specific CHARMM pair style with potentially adjusted/scale LJ parameters for 1-4 non-bonded interactions.

axel.

Thank you for the elaborating. But in AMBER, especially for aromatic rings, some 1-4 non-bonded terms of particular dihedrals are turned off. While the “special_bonds” sets a global scaling factor, is there a way to turn off some 1-4 terms in particular dihedrals of the ring system?

Andy

Thank you for the elaborating. But in AMBER, especially for aromatic rings, some 1-4 non-bonded terms of particular dihedrals are turned off. While the “special_bonds” sets a global scaling factor, is there a way to turn off some 1-4 terms in particular dihedrals of the ring system?

no. are you sure you are not confusing things here?
can you provide an exact definition of what should be excluded where and how?

axel.

Please see the following link for the AMBER prmtop instructions. https://ambermd.org/FileFormats.php
In the dihedral section, it says:

"The true atom number equals the absolute value of the number divided by three, plus one. In the case of the dihedrals, if the fourth atom is negative, this implies that the dihedral is an improper. If the third atom is negative, this implies that the end group interations are to be ignored. End group interactions are ignored, for example, in dihedrals of various ring systems (to prevent double counting of 1-4 interactions) and in multiterm dihedrals.”

Thank you for the explaining before. Now I am assured that LAMMPS won’t double count multiterm dihedrals. I am attaching a phenol molecule built from AMBER package. As you can see from the topology file, in the "%FLAG DIHEDRALS_WITHOUT_HYDROGEN“ section, you can see three dihedrals,

"0 15 -12 9 1”
"3 6 -9 12 1”
"6 9 -12 15 1”

They all have the third atom to be negative. I am changing that into readable format. Those three dihedrals are atom #1-6-5-4 ; 2-3-4-5; 3-4-5-6 with dihedral typeID 1; and, for my phenol molecule, 1 to 6 atoms are the aromatic carbon. As the AMBER description says, their 1-4 terms should be turned off. I am also attaching my mol2 file for a reference. I used GAFF for the phenol molecule.

Thank you.
Andy

top (12.1 KB)

phenol.mol2 (1.52 KB)

it is not clear from the description you pointed to what is meant by “end group interactions”. that is a very generic term.

let me reiterate the situation in LAMMPS: 1-4 non-bonded terms are guaranteed to be computed only once, since they are computed as part of the regular non-bonded interactions through walking the neighborlists. the 1-2, 1-3, 1-4 scaling factors set through special_bonds are applied. double counting of these non-bonded interactions is impossible by construction.

the exception is for the dihedral style charmm in combination with a pair style lj/charmm/coul/* where those pairwise 1-4 interactions are required to be zero in the special_bonds settings and thus not part of the pair style computation and instead computed as part of the charmm dihedral style. this feature is turned on when the 1-4 scaling factor is set to a non-zero number. the reason for this is that CHARMM may have different lennard-jones parameters for 1-4 non-bonded terms, that are stored with the pair style.
please note that this scaling factor only applies to the non-bonded interactions, NOT the dihedral itself. the scaling factor must be set to 0.5 for dihedral types in 6-rings like benzene and 0.0 for dihedral types in 5-rings or smaller, since otherwise the non-bonded interactions would be either double counted in 6-ring or included incorrectly since the same 1-4 pair of dihedral atoms is also a 1-3 pair and thus fully excluded.

since, to the best of my knowledge, amber does not have alternate LJ parameters for 1-4 interactions, they can be computed as part of the regular pair style (with the appropriate scaling factor) and then the scaling factor in the dihedral style must be 0. since you are supposed to use a non-CHARMM pair style, this is your only option anyway.

axel.