Re: [AMBER] Setting exclusions - am i overcounting ?-Last question

From: Arun Srikanth <askforarun.gmail.com>
Date: Thu, 16 Oct 2014 01:35:59 -0400

Yes thank you very much for the detailed reply. I have collected your
replies and Daniels replies on this topic, so in case I forget how things
are implemented I can look at it.

You are right. LAMMPS does calculate it correctly. I did have a look at the
code in lammps and the mailing list. So once a pair is identified, there is
no way it is double counted. Also the 1-2 weighting (0) takes precedence
over 1-3 weighting (0) and 1-3 over 1-4 (1/2) which is the same thing you
mentioned. So in the case of the ring structure the exclusion of 1-4
interaction is automatically taken care based on precedence.

This discussion was extremely helpful and has finally put my mind to peace
after a long time.

Thank you once again daniel and jason
Arun

On Thu, Oct 16, 2014 at 12:24 AM, Jason Swails <jason.swails.gmail.com>
wrote:

> On Wed, Oct 15, 2014 at 10:54 PM, Arun Srikanth <askforarun.gmail.com>
> wrote:
>
> > That was an awesome explanation for the impropers (I should have thought
> > about bonds angles taking precedence) thanks. This fact has not been
> > adressed properly anywhere atleast as far I have searched. I couldn't
> find
> > this in cornell's paper or GAFF (junmei) or may be it is obvious for
> people
> > who have experience in MD.
> >
> > I rely on parmed.py extensively for my simulations, often to check the
> > dihedral connectivity information and also for exclusions.
> >
> > Your statment :
> >
> > "However, since the
> > atom types defined in that torsion *may* be used elsewhere in your
> > molecule, you need to make sure that the nonbonded scaling factors for
> that
> > parameter is "correct".
> >
>
> ​I've thought of a clearer way of saying this. The scaling factors are
> defined in the parameter files (parm10.dat and frcmod files), with the
> default values of 1.2 and 2.0 (EEL and vdW) assigned for any torsion whose
> scaling factors are not defined.
>
> ​The only way you can know for *certain* that the end-group interactions
> for a particular torsion term will *always* be excluded from the 1-4
> nonbonded terms is if it is one of the first N-1 terms of an N-term
> torsion. Only then can you assign scaling factors of 0 in the parameter
> file. [1] The torsions you mentioned previously that was excluded (by
> virtue of being in a ring) was a single, normal torsion term defined
> between 4 atom types, so it needs to assign the "normal" scaling factors
> (normal being whatever is normal for that force field). So even if a
> torsion is only used in ring systems where those interactions are excluded,
> you don't know that for certain when you write the force field parameter
> database, so you can't set the scaling factors to 0.
>
> To take into consideration the above statement, in my simulation I have the
> > number of dihedral types equal to number of dihedrals that is each
> > dihedral has its own type (as defined by parmed.py).
>
>
> ​Note that tleap does not do this by default (nor do I think there is an
> option). There is a DIHEDRAL_PHASE, DIHEDRAL_PERIODICITY, and
> DIHEDRAL_FORCE_CONSTANT array (and SCXX_SCALE_FACTOR​
>
> ​when the scale factors are assigned)​. Each dihedral listed in
> DIHEDRALS_INC_HYDROGEN and DIHEDRALS_WITHOUT_HYDROGEN contain 5 numbers --
> the atom indexes and a pointer into the 3 (or 5) arrays I previously
> listed. In my experience, the combined lengths of the 2 dihedral arrays
> has never equalled the lengths of the 3 (or 5) parameter arrays (although
> that's not impossible).
>
> Although I am afraid
> > that the end group interactions (1-4 non bonded) of a multiterm dihedral
> > (atoms with the same type ) may be counted twice ( I do not know yet how
> > lammps calculates the 1-4 interaction for a multi term dihedral, if the
> > multiterm dihedral is defined twice as two new type of dihedrals)
> > ​
> >
>
>
> >
> > Your opinion on my above imeplementation would be helpful
> >
>
> ​​
> I
> ​ would hope that LAMMPS does it correctly. I would bet that they do,
> since the other common biomolecular force fields work the same way (like
> CHARMM, for instance). Note that using the dihedral list to build the 1-4
> nonbonded pair list (that have to be treated *specially*) is an
> implementation detail. Amber does it because it's convenient, and they
> "decorate" the dihedral lists with negatives in order to make the code
> easier to write. But you certainly don't *have* to do it this way. You
> can just look at all of the bonds, and build the 1-4 nonbonded arrays from
> unique collections of atoms that are separated by at least 3 bonds (note in
> 5-member rings, torsions are separated by 3 bonds in 1 direction and 2
> bonds in another, so they don't get added to this list).
>
> That is an equally valid, and definitely more general, way to define the
> 1-4 nonbonded pair list. Both will give you the same answer for the Amber
> force field (again, unless you used ParmEd to mess around with the dihedral
> or bond list).
>
> ​What I hope you've gathered from my tl;dr answers so far is that many of
> the file formats in Amber are intricately linked to implementation details
> in the code or relics from a code that once was (for instance, why are all
> atom indexes represented as multiples of 3 indexed from 0? why are there
> negatives in torsion terms? ...etc.) The 1-4 pair list is derived from
> the torsion list, which is why some things in the prmtop are the way they
> are. This makes some things easy (like coding the Amber force field), but
> can make other things more difficult (like supporting coarse-grained force
> fields and even the CHARMM force field). Programs that were designed to
> support multiple force fields and file types (like NAMD, OpenMM, maybe
> LAMMPS ??) may make very different implementation choices than Amber,
> because they were guided by different design principles and different
> goals. [2]
>
> ​HTH,
> Jason
>
> [1] This is actually done in tleap, but based on the same principle.
>
> [2] I've found that reading Amber code helps decode certain file formats --
> see http://ambermd.org/doc/OFF_file_format.txt for example -- and vice
> versa (e.g., I knew what to look for in the code based on how a file format
> was laid out, like with the 1-4 scaling terms defined in the prmtop). In
> many cases these file formats and some old code relics help tell the story
> of how the Amber codes have developed.
>
> --
> Jason M. Swails
> BioMaPS,
> Rutgers University
> Postdoctoral Researcher
> _______________________________________________
> 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 Oct 15 2014 - 23:00:02 PDT
Custom Search