[gmx-users] Important: Bugs in NEMD calculation
David van der Spoel
spoel at xray.bmc.uu.se
Wed Jan 19 19:13:12 CET 2011
On 2011-01-19 18.36, Xiaohu Li wrote:
> Hi, All,
> I've found a bug in the NEMD calculation for viscosity. This has
> been reviewed in /*Hess's paper at JCP 116 209 2002.*/
> The version of gromacs I'm using is the development version.
> Notice that this version correct a previous but of 4.5.3, where you uses
> NEMD, both the term V(eq. 21 on Hess's paper) and 1/eta(shear viscosity
> inverse) are supposed to be written to the *.edr file, however,
> the 4.5.3 versions does not have this. This version can be retrieved at
> *However, even this version is buggy since eta=A*rho/(V*k**2)(eq. 20
> Hess's paper)*. *I have performed simulation and has found out that the
> V and eta which are written in *edr file does not match according to the
> formula, a little further check on the source code mdebin.c under the
> directory src/mdlib shows that it is actually calculating
> eta=A*Volume/(V*k**2) where density of rho should have been used. (This
> is at line 755 of mdebin.c ).
> I hope everyone who is using this can be aware of this, if you ever
> used this code to produce data, the V is correct from *edr, however, you
> need to manuelly get your eta using the above formula.
> For the GMX developers, I hope anyone of you can correct this bug.
Thanks for pointing that out. There also is a small issue in that the
volume is computed for a rectangular box
vol = box[XX][XX]*box[YY][YY]*box[ZZ][ZZ];
dens = (tmass*AMU)/(vol*NANO*NANO*NANO);
which would be incorrect for a non-rectangular box. You should however
use a rectangular box for this kind of calculations, although this is
not enforced by grompp.
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone: +46184714205.
spoel at xray.bmc.uu.se http://folding.bmc.uu.se
More information about the gromacs.org_gmx-users