<html>
<head>
<style><!--
.hmmessage P
{
margin:0px;
padding:0px
}
body.hmmessage
{
font-size: 10pt;
font-family:Tahoma
}
--></style>
</head>
<body class='hmmessage'>
<br><br>&gt; Date: Wed, 19 Jan 2011 20:52:13 +0100<br>&gt; From: spoel@xray.bmc.uu.se<br>&gt; To: gmx-users@gromacs.org<br>&gt; Subject: Re: [gmx-users] Important: Bugs in NEMD calculation<br>&gt; <br>&gt; On 2011-01-19 20.43, Berk Hess wrote:<br>&gt; &gt;<br>&gt; &gt;<br>&gt; &gt;  &gt; Date: Wed, 19 Jan 2011 19:13:12 +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] Important: Bugs in NEMD calculation<br>&gt; &gt;  &gt;<br>&gt; &gt;  &gt; On 2011-01-19 18.36, Xiaohu Li wrote:<br>&gt; &gt;  &gt; &gt; Hi, All,<br>&gt; &gt;  &gt; &gt; I've found a bug in the NEMD calculation for viscosity. This has<br>&gt; &gt;  &gt; &gt; been reviewed in /*Hess's paper at JCP 116 209 2002.*/<br>&gt; &gt;  &gt; &gt; The version of gromacs I'm using is the development version.<br>&gt; &gt;  &gt; &gt; Notice that this version correct a previous but of 4.5.3, where you<br>&gt; &gt; uses<br>&gt; &gt;  &gt; &gt; NEMD, both the term V(eq. 21 on Hess's paper) and 1/eta(shear viscosity<br>&gt; &gt;  &gt; &gt; inverse) are supposed to be written to the *.edr file, however,<br>&gt; &gt;  &gt; &gt; the 4.5.3 versions does not have this. This version can be retrieved at<br>&gt; &gt;  &gt; &gt;<br>&gt; &gt; http://repo.or.cz/w/gromacs.git/commit/c83de86d65ce7135be6cef15e9100d7516e6d9a7<br>&gt; &gt;  &gt; &gt; *However, even this version is buggy since eta=A*rho/(V*k**2)(eq. 20<br>&gt; &gt;  &gt; &gt; Hess's paper)*. *I have performed simulation and has found out that the<br>&gt; &gt;  &gt; &gt; V and eta which are written in *edr file does not match according<br>&gt; &gt; to the<br>&gt; &gt;  &gt; &gt; formula, a little further check on the source code mdebin.c under the<br>&gt; &gt;  &gt; &gt; directory src/mdlib shows that it is actually calculating<br>&gt; &gt;  &gt; &gt; eta=A*Volume/(V*k**2) where density of rho should have been used. (This<br>&gt; &gt;  &gt; &gt; is at line 755 of mdebin.c ).<br>&gt; &gt;  &gt; &gt; I hope everyone who is using this can be aware of this, if you ever<br>&gt; &gt;  &gt; &gt; used this code to produce data, the V is correct from *edr,<br>&gt; &gt; however, you<br>&gt; &gt;  &gt; &gt; need to manuelly get your eta using the above formula.<br>&gt; &gt;  &gt; &gt; For the GMX developers, I hope anyone of you can correct this bug.<br>&gt; &gt;  &gt; &gt;<br>&gt; &gt;<br>&gt; &gt; I fixed this bug recently is the git release-4-5-patches branch.<br>&gt; &gt; You can get the fix from git and it will be in the 4.5.4 release<br>&gt; &gt; (no date yet).<br>&gt; &gt;<br>&gt; &gt;  &gt; Thanks for pointing that out. There also is a small issue in that the<br>&gt; &gt;  &gt; volume is computed for a rectangular box<br>&gt; &gt;  &gt; vol = box[XX][XX]*box[YY][YY]*box[ZZ][ZZ];<br>&gt; &gt;  &gt; dens = (tmass*AMU)/(vol*NANO*NANO*NANO);<br>&gt; &gt;  &gt;<br>&gt; &gt;  &gt; which would be incorrect for a non-rectangular box. You should however<br>&gt; &gt;  &gt; use a rectangular box for this kind of calculations, although this is<br>&gt; &gt;  &gt; not enforced by grompp.<br>&gt; &gt;<br>&gt; &gt; No.<br>&gt; &gt; That formula is correct for any triclinic box!<br>&gt; <br>&gt; But the volume of a triclinic box is the determinant of the box (maybe <br>&gt; less than the product of the XX YY ZZ).<br><br>Gromacs only allows boxes with the upper (or lower depending on how you<br>order the elements) triangle zero. This means the determinant is equal<br>to the product of the diagonal elements. This product is (slighlty) more accurate<br>and computationally less expensive than calling the det function.<br><br>Berk<br><br>&gt; <br>&gt; &gt;<br>&gt; &gt; Berk<br>&gt; &gt;<br>&gt; &gt;  &gt;<br>&gt; &gt;  &gt; &gt;<br>&gt; &gt;  &gt; &gt;<br>&gt; &gt;  &gt; &gt; Xiaohu<br>&gt; &gt;  &gt; &gt; *<br>&gt; &gt;  &gt; &gt;<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; Dept. of Cell &amp; Molec. Biol., Uppsala University.<br>&gt; &gt;  &gt; Box 596, 75124 Uppsala, Sweden. Phone: +46184714205.<br>&gt; &gt;  &gt; spoel@xray.bmc.uu.se 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://lists.gromacs.org/mailman/listinfo/gmx-users<br>&gt; &gt;  &gt; Please search the archive at<br>&gt; &gt; http://www.gromacs.org/Support/Mailing_Lists/Search before 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/Support/Mailing_Lists<br>&gt; &gt;<br>&gt; <br>&gt; <br>&gt; -- <br>&gt; David van der Spoel, Ph.D., Professor of Biology<br>&gt; Dept. of Cell &amp; Molec. Biol., Uppsala University.<br>&gt; Box 596, 75124 Uppsala, Sweden. Phone:        +46184714205.<br>&gt; spoel@xray.bmc.uu.se    http://folding.bmc.uu.se<br>&gt; -- <br>&gt; gmx-users mailing list    gmx-users@gromacs.org<br>&gt; http://lists.gromacs.org/mailman/listinfo/gmx-users<br>&gt; Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/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/Support/Mailing_Lists<br>                                               </body>
</html>