[gmx-users] rigid tetrahedral molecule
Mark Abraham
Mark.Abraham at anu.edu.au
Thu Apr 28 12:35:44 CEST 2011
On 4/28/2011 7:07 PM, Sanku M wrote:
> Thanks Mark for your comment. But, as far as degrees of freedom is
> concerned, if there is a tetrahedral molecule as I shown below with
> atom #1 ( B) being center of the tetrahedron , if we had defined the
> molecules in terms of bonds and angles ( in stead of constraints ), we
> would have 4 bonds ( 1-2,1-3,1-4,1-5 ) and 6 angles (
> 2-1-3,2-1-4,2-1-5,3-1-4,3-1-5 , 4-1-5 )at 109.45 degrees , thus total
> ( 4+6)=10 descriptors .
Different things happen in each case. A redundant but consistent bonded
interaction adds (at most) a small contribution to the energies and
forces. Consider what the effect of a suitably-constructed H-H-C angle
function of a methylene in an equilibrium configuration would be...
nothing significant.
However, an implementation of coupled constraints so that they are all
mutually satisifed in few iterations is hard enough without having to
cope with redundancy also.
> In that case, the simulation goes fine . But, then wonder, why, for
> constraints, we can put only 9 descriptors .
> Also, I have tried all possible combination to get 9 constraints, but
> , in each case, gromacs crashes. However, defining 4 bonds and 6
> angles using high force-constant lets the simulation go fine. So, on
> a second thought, I wonder whether defining 9 bond-constraints to get
> a stable tetrahedral molecule is at all possible.
So you need to read the LINCS and/or SHAKE literature :-) You say these
other codes can do rigid-body simulations, so there must be literature
on how to do it.
> Do I need to go for some sort of compromise, where 4 bonds
> ( 1-2,1-3,1-4,1-5 ) are defined as constraints and 6 angles with high
> force-constant are defined to maintain tetrahedral nature ?
Maybe.
> I am not sure what should be the best possible 9 constraints .
Don't know. I would start trying four B-F bonds and five F-F bonds, but
it shouldn't matter.
Mark
> 1 opls_1014 1 BF4 B 1 0.8276
>>> 2 opls_1015 1 BF4 F1 1 -0.4569
>>> 3 opls_1015 1 BF4 F2 1 -0.4569
>>> 4 opls_1015 1 BF4 F3 1 -0.4569
>>> 5 opls_1015 1 BF4 F4 1 -0.4569
>
> ------------------------------------------------------------------------
> *From:* Mark Abraham <Mark.Abraham at anu.edu.au>
> *To:* Discussion list for GROMACS users <gmx-users at gromacs.org>
> *Sent:* Thu, April 28, 2011 1:00:42 AM
> *Subject:* Re: [gmx-users] rigid tetrahedral molecule
>
> On 4/28/2011 3:49 PM, Sanku M wrote:
>> I went through the LINCS manual . But, I am still struggling with
>> coming up with the idea of putting correct constraint to maintain the
>> rigidity of tetrahedral molecule . I seem to understand from your
>> suggestion that the tetrahedral can be seen as a combination of 4
>> coupled triangles.( or am I still wrong about it ?)
>
> I think it's seven coupled triangles, but that's not a relevant way to
> think about it. You have 5 atoms, so 3N-6 means 9 degrees of freedom,
> so you need 9 independent descriptors of relative atomic positions, so
> 9 constraints. You were trying 10 and 5.
>
> Because those constraints form triangles, you may need to take care
> with LINCS to get a stable simulation. Read up on the details here, I
> don't remember them.
>
>> In that case, am I supposed to use multiple settle to keep the
>> molecule in a tetrahedral fashion ? I am sorry but if you can explain
>> it in bit more details, I might get the point.
>
> Forget about SETTLE - it was just an example to illustrate that this
> is not easy to do right. There's a specialised algorithm for rigid
> water (with three coupled bond constraints) because it is fiddly to
> get such things right (and fast).
>
> Mark
>
>> ------------------------------------------------------------------------
>> *From:* Mark Abraham <Mark.Abraham at anu.edu.au>
>> *To:* Discussion list for GROMACS users <gmx-users at gromacs.org>
>> *Sent:* Wed, April 27, 2011 11:14:34 PM
>> *Subject:* Re: [gmx-users] rigid tetrahedral molecule
>>
>> On 4/28/2011 1:54 PM, Sanku M wrote:
>>> Hi,
>>> I tried to keep the geometry of the BF4 fixed by using constraints
>>> using lincs. But , unfortunately, my simulation is crashing
>>> immediately and if I try minimization with only 2 molecules, it
>>> provides a lot of LINCS warning and generate a lot of step*.pdb file
>>> . If I try to visualize the minimized snapshot in VMD, it looks like
>>> all the distances I tried to constrain decreased drastically.
>>> Finally, trying MD run with this "minimized" configuration results
>>> in crashing due to bad contacts.
>>>
>>> I am sure I am doing something wrong and it might be that my itp
>>> file is wrong . So any help will be highly appreciated.
>>> Here is the details of what I did.
>>>
>>> The geometry of the molecule is tetrahedral with B at the
>>> center and 4 F atoms is surrounding it in a tetrahedral manner.
>>> I first generated a itp file for BF4 which is shown below: I first
>>> got the LJ parameters and charges for B and F atom and put them in
>>> ffoplsnb.itp file as new atom types opls_1014 and opls_1015
>>> . Initially I tried to put contsraint along all bonds ( i.e among F
>>> atoms as well ). But, grompp provides warning that number of
>>> constraint is more than number of degrees of freedom. So, I reduced
>>> number of constraints by only putting constraint among B and F. But,
>>> it did not work either.
>>
>> Sure, you need as many constraints as
>> http://en.wikipedia.org/wiki/Degrees_of_freedom_%28physics_and_chemistry%29
>>
>> You should also do your homework about using LINCS and coupled
>> triangles of constraints, as I suggested last time.
>>
>> Mark
>>
>>> Here is the .itp file I wrote for rigid BF4 . It will be great if
>>> someone can point me what I am doing wrong.
>>>
>>> [ moleculetype ]
>>> ; molname nrexcl
>>> BF4 3
>>>
>>> [ atoms ]
>>> #ifdef _FF_OPLS
>>> 1 opls_1014 1 BF4 B 1 0.8276
>>> 2 opls_1015 1 BF4 F1 1 -0.4569
>>> 3 opls_1015 1 BF4 F2 1 -0.4569
>>> 4 opls_1015 1 BF4 F3 1 -0.4569
>>> 5 opls_1015 1 BF4 F4 1 -0.4569
>>> #endif
>>> [ constraints ]
>>> 1 2 1 0.146
>>> 1 3 1 0.146
>>> 1 4 1 0.146
>>> 1 5 1 0.146
>>> ; 2 3 1 0.238
>>> ; 2 4 1 0.238
>>> ; 2 5 1 0.238
>>> ; 3 4 1 0.238
>>> ; 3 5 1 0.238
>>> ; 4 5 1 0.238
>>>
>>> [ exclusions ]
>>> 1 2 3 4 5
>>> 2 1 3 4 5
>>> 3 1 2 4 5
>>> 4 1 2 3 5
>>> 5 1 2 3 4
>>>
>>>
>>> ------------------------------------------------------------------------
>>> *From:* Mark Abraham <Mark.Abraham at anu.edu.au>
>>> *To:* Discussion list for GROMACS users <gmx-users at gromacs.org>
>>> *Sent:* Wed, April 27, 2011 8:39:23 PM
>>> *Subject:* Re: [gmx-users] rigid tetrahedral molecule
>>>
>>> On 4/28/2011 11:25 AM, Justin A. Lemkul wrote:
>>> >
>>> >
>>> > Sanku M wrote:
>>> >> Hi,
>>> >> I am interested in simulating a anionic molecule BF4(-) ( Boron
>>> tetrafluoride). In the paper which developed the parameters for
>>> this molecule, it is mentioned that it has been used as 'rigid'
>>> molecule i.e the molecule only has non-bonding interaction but there
>>> was no intramolecular motion as the geometry was fixed.
>>> >> I am trying to simulate this molecule in gromacs treating it as
>>> rigid. But, I was looking for best way to 'rigidify' this molecule.
>>> >>
>>> >> I was wondering whether using LINCS to constrain all B-F and F-F
>>> bonds will be good enough . Or, Should I use virtual sites ? If I
>>> really need to use virtual site, will it be something like TIP5P
>>> water model ?
>>> >> Can someone suggest the best wayout ?
>>> >>
>>> >
>>> > Constraints should do the trick, but probably the best approach is
>>> to simply contact the authors who developed the model and ask how
>>> they did it. Then you know you're exactly reproducing what they did.
>>>
>>> Yep.
>>>
>>> Be aware that the coupled constraints make life tricky, and you
>>> should read up in the manual and literature for how best to use
>>> P-LINCS in such cases. Algorithms like SETTLE for rigid water exist
>>> for a reason...
>>>
>>> Mark
>>>
>>>
>>>
>>> -- gmx-users mailing list gmx-users at gromacs.org
>>> <mailto:gmx-users at gromacs.org>
>>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>>> Please search the archive at
>>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>>> Please don't post (un)subscribe requests to the list. Use the www
>>> interface or send it to gmx-users-request at gromacs.org
>>> <mailto:gmx-users-request at gromacs.org>.
>>> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20110428/fa5239c4/attachment.html>
More information about the gromacs.org_gmx-users
mailing list