Re: [gmx-users] SHAKE for water
On 10/30/13 12:00 PM, Guillaume Chevrot wrote: Hi, I was wondering if there is an option or a way to use the algorithm SHAKE for the water molecules (instead of SETTLE). Never tried it, but I would assume define = -DFLEXIBLE constraint-algorithm = shake would do the trick. Of course, that assumes that your water model was not parametrized to stay rigid with SETTLE (which most were) and that you're OK using SHAKE for everything else in your system, since you can only use one constraint algorithm. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- gmx-users mailing listgmx-users@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-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] SHAKE for water
Hi, 2013/10/30 Justin Lemkul jalem...@vt.edu On 10/30/13 12:00 PM, Guillaume Chevrot wrote: Hi, I was wondering if there is an option or a way to use the algorithm SHAKE for the water molecules (instead of SETTLE). Never tried it, but I would assume define = -DFLEXIBLE constraint-algorithm = shake thanks for the trick, I'll try. would do the trick. Of course, that assumes that your water model was not parametrized to stay rigid with SETTLE (which most were) Does that mean changing [ settle ] by [ constraints ] in the file spce.itp is sufficient? Thanks again! Guillaume and that you're OK using SHAKE for everything else in your system, since you can only use one constraint algorithm. -Justin -- ==** Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalemkul@outerbanks.umaryland.**edu jalem...@outerbanks.umaryland.edu | (410) 706-7441 ==** -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/**mailman/listinfo/gmx-usershttp://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/** Support/Mailing_Lists/Searchhttp://www.gromacs.org/Support/Mailing_Lists/Searchbefore posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/**Support/Mailing_Listshttp://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@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-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] SHAKE for water
On 10/30/13 1:46 PM, Guillaume Chevrot wrote: Hi, 2013/10/30 Justin Lemkul jalem...@vt.edu On 10/30/13 12:00 PM, Guillaume Chevrot wrote: Hi, I was wondering if there is an option or a way to use the algorithm SHAKE for the water molecules (instead of SETTLE). Never tried it, but I would assume define = -DFLEXIBLE constraint-algorithm = shake thanks for the trick, I'll try. would do the trick. Of course, that assumes that your water model was not parametrized to stay rigid with SETTLE (which most were) Does that mean changing [ settle ] by [ constraints ] in the file spce.itp is sufficient? Sufficient for what? What type of replacement are you making? If you're just trying to use SHAKE, you don't have to modify the topology at all; the .mdp options are sufficient to produce the behavior you want. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- gmx-users mailing listgmx-users@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-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] shake for water
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: sp...@xray.bmc.uu.se To: gmx-users@gromacs.org Subject: Re: [gmx-users] shake for water Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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@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-requ...@gromacs.org. Can't post? Read http
Re: [gmx-users] shake for water
David Mobley wrote: 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? Of course you are welcome to send a bugzilla. You can also use g_bond to plot a bond histogram of water. (g_bond -blen 0.1 -tol 1e-5 or something along those lines) Thanks, David Berk Berk Date: Thu, 5 Feb 2009 22:41:47 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@gromacs.org Subject: Re: [gmx-users] shake for water Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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@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-requ...@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. sp
Re: [gmx-users] shake for water
Berk, Another question, just to be sure. Have you actually checked that the other code really gets the distances right up to 1e-12? Yes. Berk Date: Thu, 5 Feb 2009 12:03:21 -0600 From: dmob...@gmail.com To: gmx-users@gromacs.org Subject: [gmx-users] shake for water 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. Thanks, David ___ gmx-users mailing list gmx-users@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-requ...@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 listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php ___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] shake for water
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. 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? Sorry, let me clarify. I was not that clear on the difference between settle and shake. The issue is that I am using SETTLE to constrain the water and I am only getting distances right to 1e-6 though I am running in double precision. David Berk Thanks, David ___ gmx-users mailing list gmx-users@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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 listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php ___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] shake for water
Hi, I don't agree. It uses the small 1+a approximation for the square. Also mdrun prints the rmsd determined with independent code, which is consistent with the (correct) tolerance. Er, yes, I have not turned on SHAKE -- I'm trying to use LINCS for bonds to hydrogen within my molecule and SETTLE for the waters. Sorry for the confusion in my original note. Anyway -- any insights as to why I'd only be getting distances correct to 1e-6 with SETTLE? Thanks, David Berk Date: Thu, 5 Feb 2009 22:41:47 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@gromacs.org Subject: Re: [gmx-users] shake for water Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se
RE: [gmx-users] shake for water
Date: Mon, 9 Feb 2009 12:00:52 -0600 Subject: Re: [gmx-users] shake for water From: dmob...@gmail.com To: gmx-users@gromacs.org Hi, I don't agree. It uses the small 1+a approximation for the square. Also mdrun prints the rmsd determined with independent code, which is consistent with the (correct) tolerance. Er, yes, I have not turned on SHAKE -- I'm trying to use LINCS for bonds to hydrogen within my molecule and SETTLE for the waters. Sorry for the confusion in my original note. Anyway -- any insights as to why I'd only be getting distances correct to 1e-6 with SETTLE? Thanks, David 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. Berk Berk Date: Thu, 5 Feb 2009 22:41:47 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@gromacs.org Subject: Re: [gmx-users] shake for water Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx
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. Thanks, David ___ gmx-users mailing listgmx-users@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-requ...@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. sp...@xray.bmc.uu.sesp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
RE: [gmx-users] shake for water
Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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. 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 listgmx-users@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php _ Express yourself instantly with MSN Messenger! Download today it's FREE! http://messenger.msn.click-url.com/go/onm00200471ave/direct/01/___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
RE: [gmx-users] shake for water
Hi, Another question, just to be sure. Have you actually checked that the other code really gets the distances right up to 1e-12? Berk Date: Thu, 5 Feb 2009 12:03:21 -0600 From: dmob...@gmail.com To: gmx-users@gromacs.org Subject: [gmx-users] shake for water 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. Thanks, David ___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php _ Express yourself instantly with MSN Messenger! Download today it's FREE! http://messenger.msn.click-url.com/go/onm00200471ave/direct/01/___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] shake for water
Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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 listgmx-users@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-requ...@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. sp...@xray.bmc.uu.sesp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
RE: [gmx-users] shake for water
Hi, I don't agree. It uses the small 1+a approximation for the square. Also mdrun prints the rmsd determined with independent code, which is consistent with the (correct) tolerance. Berk Date: Thu, 5 Feb 2009 22:41:47 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@gromacs.org Subject: Re: [gmx-users] shake for water Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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 listgmx-users@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo
RE: [gmx-users] shake for water
Hi David, I just checked with shake tolerance 1e-10 I really get deviations of 1e-10. First I forgot to turn on Shake and got the default Lincs, which with the default settings gave me roughly the accuracy you reported. (you did not by coincidence make the same mistake?) Berk From: g...@hotmail.com To: gmx-users@gromacs.org Subject: RE: [gmx-users] shake for water Date: Thu, 5 Feb 2009 23:02:31 +0100 Hi, I don't agree. It uses the small 1+a approximation for the square. Also mdrun prints the rmsd determined with independent code, which is consistent with the (correct) tolerance. Berk Date: Thu, 5 Feb 2009 22:41:47 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@gromacs.org Subject: Re: [gmx-users] shake for water Berk Hess wrote: Date: Thu, 5 Feb 2009 19:35:09 +0100 From: sp...@xray.bmc.uu.se To: gmx-users@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@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-requ...@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. sp...@xray.bmc.uu.se sp...@gromacs.org http://folding.bmc.uu.se ___ gmx-users mailing list gmx-users@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-requ...@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 listgmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
