[gmx-users] 1,4 interactions

[Justin A. Lemkul]

Gavin Melaugh wrote:
> Hi Justin
> 
> I just wanted to say a big thanks for your time and patience
> 

No problem.  This exchange has made me realize that what's in the manual, while 
complete, lacks a bit of clarity.  I'm going to make some updates to the 
documentation for a future release.  Hopefully this will clarify some issues.

-Justin

> Gavin
> 
> Justin A. Lemkul wrote:
>>
>> Gavin Melaugh wrote:
>>> Right, I think the penny has finally dropped. On the back of what you
>>> said.
>>>
>>> 1)In my topologies I have always listed say a [bond] directive, with
>>> bond indices,function, k, and r0. I could however just list the indices
>>> and the function, and then have the corresponding parameters in the
>>> [bond types] directive, which save me associating the parameters to each
>>> individual bond ?
>>>
>> Yes.
>>
>>> 2) In reponse to what yo said finally "If you list specific pairs of
>>> atoms in a [pairs] directive, then grompp will search for the
>>> appropriate parameters listed in a [pairstypes] directive.  If none are
>>> found, and "gen-pairs = yes," grompp generates the appropriate
>>> parameters for you." Are the parameters generated from the [atomtypes]
>>> directive?
>>>
>> Yes, in conjunction with FudgeLJ and FudgeQQ.
>>
>> -Justin
>>
>>> Gavin
>>>
>>> Justin A. Lemkul wrote:
>>>> Gavin Melaugh wrote:
>>>>> A few points then,
>>>>>
>>>>> 1) So a [pairs] directive would include all the atom indices of the
>>>>> 1,4
>>>>> interactions and the coressponding values of sigma and epsilon?,
>>>>> which I
>>>>> would do manually.
>>>>>
>>>> No.  A [pairs] directive simply lists the pairs of atoms for which 1,4
>>>> interactions should be calculated.  No parameters are listed here.
>>>> Format:
>>>>
>>>> atom_i atom_j func_type
>>>>
>>>>> 2) In the manual it says ''gen-pairs = yes : generate 1,4 interactions
>>>>> which are not present in the pair list from normal LJ parameters using
>>>>> Fudge LJ". Sorry if I am coming across stupid, but that seems
>>>>> different
>>>>> to what you're saying?
>>>>>
>>>> Not at all.  The parameters for a specific interaction are given in a
>>>> [pairtypes] directive.  If these parameters are not present and
>>>> "gen-pairs = yes" then they are calculated.  If parameters for a
>>>> specific pair are not present, and "gen-pairs = no" then grompp
>>>> generates a fatal error.
>>>>
>>>>> 3) Also am Icorrect about the normal nb-LJ coming from the [atom types
>>>>> directive]
>>>>>
>>>> Yes.
>>>>
>>>>> 4) I suppose a very obvious question would be; What is the difference
>>>>> between a [pairtypes] directive and [pairs] directive (bonds,
>>>>> bondtypes etc)
>>>>>
>>>> A [*type] directive lists the parameters that should be a applied to
>>>> whatever interaction type might be found.  For instance, a [bondtypes]
>>>> directive contains default parameters for all the bonds you might
>>>> find, like C-C.  If you then have a [bonds] directive that lists:
>>>>
>>>> ; ai  aj  func
>>>>    1   2    1
>>>>
>>>> ...and atoms 1 and 2 are both C, then grompp finds the corresponding
>>>> parameters in the [bondtypes] directive and applies them to that bond,
>>>> a process which it carries out for all the rest of the interaction
>>>> types, including pairs.  If you list specific pairs of atoms in a
>>>> [pairs] directive, then grompp will search for the appropriate
>>>> parameters listed in a [pairstypes] directive.  If none are found, and
>>>> "gen-pairs = yes," grompp generates the appropriate parameters for you.
>>>>
>>>> -Justin
>>>>
>>>>> Gavin
>>>>>
>>>>> P.S I appreciate you using your time to explain this
>>>>> Justin A. Lemkul wrote:
>>>>>> Gavin Melaugh wrote:
>>>>>>> Hi Justin
>>>>>>>
>>>>>>> Thanks for replying. I have checked the energy file and there are
>>>>>>> no 1,4
>>>>>>> terms, which is what I expected. I should clarify that my
>>>>>>> parameters are
>>>>>>> taken from the OPLS ua forcefield but I have relabelled them for my
>>>>>>> own
>>>>>>> model, which is not a common molecule. Therefore all the info on the
>>>>>>> molecule is in .top and obviously the .gro file. i.e. there are no
>>>>>>> #include statements. Now I take it that the non-bonded LJ are taken
>>>>>>> from
>>>>>>> the [atom types] directive? The question is then how and where do I
>>>>>>> state the 1,4 terms. I always thought that when gen-pairs = yes,
>>>>>>> then
>>>>>>> gromacs would generate a list of all the 1,4 interaction types?
>>>>>>>
>>>>>> You don't explicitly state 1-4 terms with OPLS-AA.  While you don't
>>>>>> need a [pairtypes] directive, you do need a [pairs] directive to
>>>>>> dictate the 1-4 interactions.  Setting "gen-pairs = yes" generates
>>>>>> the
>>>>>> parameters, such that a [pairtypes] directive is unnecessary.
>>>>>>
>>>>>> -Justin
>>>>>>
>>>>>>> Gavin
>>>>>>>
>>>>>>> Justin A. Lemkul wrote:
>>>>>>>> Gavin Melaugh wrote:
>>>>>>>>> Hi all
>>>>>>>>>
>>>>>>>>> Can someone please get back to me on this. I have generated my own
>>>>>>>>> topology file from OPLS parameters by hand with the following
>>>>>>>>> format. I
>>>>>>>>> am now concerned that I haven't been calculating 1,4 interactions,
>>>>>>>>> which
>>>>>>>>> I thought were generated by setting gen-pairs = yes. I have read
>>>>>>>>> section
>>>>>>>>> 5 of the manual several times but am still unclear about this.
>>>>>>>>>
>>>>>>>> You can answer this question by analyzing your .edr file: are there
>>>>>>>> LJ-14 and/or Coul-14 terms?  If so, then you're calculating 1-4
>>>>>>>> interactions.  Seems to me that, for OPLS, your approach is
>>>>>>>> fine.  It
>>>>>>>> wasn't clear from your earlier posts that you were using OPLS,
>>>>>>>> hence
>>>>>>>> the inability to know if what you were doing was right or not.
>>>>>>>> Different force fields have different underlying mechanics, so the
>>>>>>>> use
>>>>>>>> of [nonbond_params] and/or [pairs]/[pairtypes] is somewhat
>>>>>>>> dependent
>>>>>>>> upon your model.
>>>>>>>>
>>>>>>>> -Justin
>>>>>>>>
>>>>>>>>> Cheers
>>>>>>>>>
>>>>>>>>> Gavin
>>>>>>>>>
>>>>>>>>> parameter level
>>>>>>>>>> [defaults]
>>>>>>>>>> ; nbfunc    comb-rule     gen-pairs        fudgeLJ     fudgeQQ
>>>>>>>>>>      1         3              yes            0.5         0.5
>>>>>>>>>>
>>>>>>>>>> [atomtypes]
>>>>>>>>>> ;type     mass           charge      ptype     sigma(nm) 
>>>>>>>>>> epsilon(kjmol-1)
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> ;Molecular level
>>>>>>>>>> [moleculetype]
>>>>>>>>>> ;       name         nrexcl
>>>>>>>>>>         tertpentane      3
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> [atoms]
>>>>>>>>>> ; atomnr  type   resnr  residue   name    cgnr     charge    
>>>>>>>>>> mass
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> [bonds]
>>>>>>>>>> ;  ai    aj   func   b0(nm)       kb(kjmol-1 nm-2)
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> [angles]
>>>>>>>>>> ;  ai    aj    ak   func  theta(deg)    ktheta(kjmol-1 rad-2)
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> [dihedrals]
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> ;System level
>>>>>>>>>> [system]
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> [molecules]
>>>>>>>>>> ; moleculename     number
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>   
>>>
> 
> 

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================