[gmx-users] Simulating Multiple Solute Particles

[Justin Lemkul]

On 11/10/14 3:07 AM, Nathan K Houtz wrote:
> Thank you again, Dr. Lemkul! I went with the first option because my advisor wants me to use tip3p water. Unfortunately, I immediately got myself stuck on the next step, but this time I think the solution may not be as trivial. The command to output em.tpr seems to work without any errors. The command,
>

Well, strictly speaking, either option works because the parameters of TIP3P are 
not force field-dependent, just translated between different conventions for 
nonbonded parameters, but that's just a pedantic issue.

>> gmx mdrun -deffnm em
>
> gives me a lot of errors. (It stops after 1000). Here's a small sample:
>
>> Steepest Descents:
>>    Tolerance (Fmax)   =  1.00000e+03
>>    Number of steps    =        50000
>>
>> Step 1, time 0.001 (ps)  LINCS WARNING
>> relative constraint deviation after LINCS:
>> rms 0.373893, max 0.619462 (between atoms 1747 and 1748)
>> bonds that rotated more than 30 degrees:
>> atom 1 atom 2  angle  previous, current, constraint length
>>       4      5   74.8    0.1130   0.0501      0.1204
>>      11     12   75.0    0.1132   0.0510      0.1204
>>      18     19   74.9    0.1129   0.0497      0.1204
>>      25     26   75.6    0.1126   0.0488      0.1204
>>      32     33   74.9    0.1132   0.0513      0.1204
>>      39     40   75.7    0.1129   0.0503      0.1204
>
> etc. There are 7 atoms (the methyl group is a united atom, otherwise there would be 10) in each TTA molecule, so you'll notice that it's complaining about the bond between the 4th and 5th atoms on every one of them. My best guess is that my topology file has insufficient or incorrect constraints, allowing this particular bond to rotate. The entire molecule is supposed to be rigid. I did constrain all the distances, but not the angles or dihedrals. One of the example topology files I used did the same: constrained distances but not angles or dihedrals. I was skeptical because it seems to allow some degrees of freedom but I assumed gromacs somehow took care of it. Even in this instance, it doesn't make sense to me that only one of the bonds is problematic if none of them are constrained. So, I'm unsure of what exactly is being left out. The bond in question is between the double bonded C=O. My molecule looks like this:
>
> H              O
>    \          //
> H - C ≡ C - C
>    /          \
> H             O-H
>
> My topology file is available from a previous message in this thread (minus the last force field correction), linked here from the archive: https://www.mail-archive.com/gromacs.org_gmx-users@maillist.sys.kth.se/msg07776.html
>
> If anyone knows directly how to fix the error, I would greatly appreciate it, but helping me understand the error would also be very useful: Does "relative constraint deviation" mean that a constraint has been broken, or that there is no constraint where there should be one? Also, when it's saying the bonds are rotating more than 30 degrees, is that the bond angle changing or does it mean that a dihedral is changing by that much?
>

Your system is http://www.gromacs.org/Documentation/Terminology/Blowing_Up 
because the topology is probably unstable.  Triple bonds need to be handled 
using virtual sites.  See previous discussions about linear construction 
(examples are acetonitrile and others) and my tutorial on doing just this.

-Justin

> Thanks again in advance - I really appreciate the help!
> N.H.
>
> ----- Original Message -----
> From: "Justin Lemkul" <jalemkul at vt.edu>
> To: gmx-users at gromacs.org
> Sent: Tuesday, November 4, 2014 7:45:39 AM
> Subject: Re: [gmx-users] Simulating Multiple Solute Particles
>
>
>
> On 11/4/14 12:45 AM, Nathan K Houtz wrote:
>> Thanks for your reply. However, I'm still confused. I thought that the command:
>>
>> #include "oplsaa.ff/tip3p.itp"
>>
>> is a call to the opls-aa force field. If this is not the correct way to include the force field parameters, how should I do that?
>>
>
> That line does not call a force field, it calls a topology for TIP3P, which
> makes use of OPLS-AA parameters (well, more correctly TIP3P parameters as
> translated into OPLS-AA atom types and functional form).  Your options are
> either (1) #include "oplsaa.ff/forcefield.itp" at the start of the topology and
> remove the [defaults] directive in your TETR topology or (2) translate TIP3P
> into atom types that are self-contained in your topology.
>
> -Justin
>

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

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

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