<html>
<head>
<style>
.hmmessage P
{
margin:0px;
padding:0px
}
body.hmmessage
{
font-size: 10pt;
font-family:Verdana
}
</style>
</head>
<body class='hmmessage'>
Hi David,<br><br>I just checked with shake tolerance 1e-10 I really get deviations of 1e-10.<br>First I forgot to turn on Shake and got the default Lincs,<br>which with the default settings gave me roughly the accuracy you reported.<br>(you did not by coincidence make the same mistake?)<br><br>Berk<br><br><hr id="stopSpelling">From: gmx3@hotmail.com<br>To: gmx-users@gromacs.org<br>Subject: RE: [gmx-users] shake for water<br>Date: Thu, 5 Feb 2009 23:02:31 +0100<br><br>



<style>
.ExternalClass .EC_hmmessage P
{padding:0px;}
.ExternalClass body.EC_hmmessage
{font-size:10pt;font-family:Verdana;}
</style>


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