[gmx-users] shake for water

[David Mobley]

Berk,
> Ah, now we are talking about something completely different.
> Settle gets me at least 1e-10.
> Are you sure your problem is in the water and not in the molecule with
> LINCS?
> For LINCS you will need to change lincs_order to at least 8 and lincs_iter
> to something
> probably something around 8.

I am using trjconv with -ndec to output the coordinates to gro files
to 12 decimal places, then computing the distances between the atoms
in the water from these coordinates. So yes, I am sure I am looking at
the water molecules. The distances in question are only right to
~1e-6.

Should I send input files?

Thanks,
David


> Berk
>
>
>>
>> >
>> > Berk
>> >
>> >> Date: Thu, 5 Feb 2009 22:41:47 +0100
>> >> From: spoel at xray.bmc.uu.se
>> >> To: gmx-users at gromacs.org
>> >> Subject: Re: [gmx-users] shake for water
>> >>
>> >> Berk Hess wrote:
>> >> >
>> >> >
>> >> > > Date: Thu, 5 Feb 2009 19:35:09 +0100
>> >> > > From: spoel at xray.bmc.uu.se
>> >> > > To: gmx-users at gromacs.org
>> >> > > Subject: Re: [gmx-users] shake for water
>> >> > >
>> >> > > David Mobley wrote:
>> >> > > > All,
>> >> > > >
>> >> > > > A quick question on constraints... I'm using TIP4P-Ew in gromacs
>> >> > > > 3.3.2
>> >> > > > and am concerned with reproducing energies from another code very
>> >> > > > precisely for several specific snapshots. I am doing a zero-step
>> >> > > > mdrun
>> >> > > > of a setup with one small molecule and two tip4p-ew water
>> >> > > > molecules.
>> >> > > >
>> >> > > > Anyway, I have set the shake tolerance to 1e-12 in the mdp file,
>> >> > > > but
>> >> > > > to my surprise the internal water distances are good only to
>> >> > > > 1e-06
>> >> > > > and
>> >> > > > 1e-07. Is this expected behavior? Note that I am running in
>> >> > > > double
>> >> > > > precision. I assumed that, er, the distances should converge to
>> >> > > > the
>> >> > > > shake tolerance.
>> >> > > Well, the documentation might be lacking, but the code tells the
>> >> > > truth.
>> >> > > It seems that the tolerance is used on the distance squared, which
>> >> > > is
>> >> > > consistent with your observation of a precision of 1e-6. So try
>> >> > > 1e-24.
>> >> >
>> >> > No, it is not the square.
>> >> > The code does a small 'a' approximation: (1+a)^2=1+2a+a^2 is approx
>> >> > 1+2a.
>> >> > I have also tested/benchmark shake in gromacs for my first lincs
>> >> > paper
>> >> > and the plincs paper and it always behaved the way I thought it
>> >> > would.
>> >>
>> >>
>> >> First we compute the inverse square of the shake distance dA in tt[ll]
>> >> for each shake pair:
>> >>
>> >> if (bFEP)
>> >> toler = sqr(L1*ip[type].shake.dA + lambda*ip[type].shake.dB);
>> >> else
>> >> toler = sqr(ip[type].shake.dA);
>> >> dist2[ll] = toler;
>> >> tt[ll] = 1.0/(toler*tol2);
>> >> }
>> >>
>> >> Then in the shake iteration we compute the difference between the
>> >> squared distances (variable diff below):
>> >>
>> >> tx = xp[ix]-xp[jx];
>> >> ty = xp[iy]-xp[jy];
>> >> tz = xp[iz]-xp[jz];
>> >> rpij2 = tx*tx+ty*ty+tz*tz;
>> >> toler = dist2[ll];
>> >> diff = toler-rpij2;
>> >>
>> >> Now we multiply the diff with tt[ll], in other words we get
>> >> iconv = (1-rpij2/dist2)/(2 tol)
>> >>
>> >> /* iconv is zero when the error is smaller than a bound */
>> >> iconv = fabs(diff)*tt[ll];
>> >>
>> >> In other words, the tolerance operates on the squared distance.
>> >>
>> >>
>> >>
>> >>
>> >> >
>> >> > The problem is not that you need more iterations than 1000?
>> >> >
>> >> > And why not use settle that does it right at full precision at once?
>> >> >
>> >> > Berk
>> >> >
>> >> >
>> >> > >
>> >> > > >
>> >> > > > Thanks,
>> >> > > > David
>> >> > > > _______________________________________________
>> >> > > > gmx-users mailing list gmx-users at gromacs.org
>> >> > > > http://www.gromacs.org/mailman/listinfo/gmx-users
>> >> > > > Please search the archive at http://www.gromacs.org/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.
>> >> > > > Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>> >> > >
>> >> > >
>> >> > > --
>> >> > > David van der Spoel, Ph.D., Professor of Biology
>> >> > > Molec. Biophys. group, Dept. of Cell & Molec. Biol., Uppsala
>> >> > > University.
>> >> > > Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. Fax:
>> >> > > +4618511755.
>> >> > > spoel at xray.bmc.uu.se spoel at gromacs.org http://folding.bmc.uu.se
>> >> > > _______________________________________________
>> >> > > gmx-users mailing list gmx-users at gromacs.org
>> >> > > http://www.gromacs.org/mailman/listinfo/gmx-users
>> >> > > Please search the archive at http://www.gromacs.org/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.
>> >> > > Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>> >> >
>> >> >
>> >> > ------------------------------------------------------------------------
>> >> > Express yourself instantly with MSN Messenger! MSN Messenger
>> >> > <http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/>
>> >> >
>> >> >
>> >> >
>> >> > ------------------------------------------------------------------------
>> >> >
>> >> > _______________________________________________
>> >> > gmx-users mailing list gmx-users at gromacs.org
>> >> > http://www.gromacs.org/mailman/listinfo/gmx-users
>> >> > Please search the archive at http://www.gromacs.org/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.
>> >> > Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>> >>
>> >>
>> >> --
>> >> David van der Spoel, Ph.D., Professor of Biology
>> >> Molec. Biophys. group, Dept. of Cell & Molec. Biol., Uppsala
>> >> University.
>> >> Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. Fax: +4618511755.
>> >> spoel at xray.bmc.uu.se spoel at gromacs.org http://folding.bmc.uu.se
>> >> _______________________________________________
>> >> gmx-users mailing list gmx-users at gromacs.org
>> >> http://www.gromacs.org/mailman/listinfo/gmx-users
>> >> Please search the archive at http://www.gromacs.org/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.
>> >> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>> >
>> > ________________________________
>> > Express yourself instantly with MSN Messenger! MSN Messenger
>> > ________________________________
>> > See all the ways you can stay connected to friends and family
>> > _______________________________________________
>> > gmx-users mailing list gmx-users at gromacs.org
>> > http://www.gromacs.org/mailman/listinfo/gmx-users
>> > Please search the archive at http://www.gromacs.org/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.
>> > Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>> >
>> _______________________________________________
>> gmx-users mailing list gmx-users at gromacs.org
>> http://www.gromacs.org/mailman/listinfo/gmx-users
>> Please search the archive at http://www.gromacs.org/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.
>> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>
> ________________________________
> Express yourself instantly with MSN Messenger! MSN Messenger
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/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.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>