[gmx-users] Virial from non-bonded forces
David van der Spoel
spoel at xray.bmc.uu.se
Tue Feb 10 18:52:06 CET 2009
Omololu Akin-Ojo wrote:
> I have two main issues:
> 1. In Appendix B of the GROMACS user manual, there is a discussion of
> the calculation of the virial. The top of the second page of the
> appendix (just before B.1.2) states, "In a triclinic system there are 27
> possible images of i, ...".
> My question is: In this case, there should also be 27 possible forces on
> particle j due to i. Granted, only the force corresponding to the
> shortest i-j distance will be non-zero but when the forces
> are multiplied by \delta_i, (in Eq. B.11) the correct force (out of all
> the 27 possible ones) should be multiplied by the corresponding \delta_i.
> In other words, if we label each of the 27 possible \delta_i by n
> ((\delta_i)^n), we should do the same with the forces F_i ((F_i)^n) and
> the sum in B.11 should also contain a sum over n.
> For ver. 3.3.1, I see that in src/mdlib/calcvir.c each of the different
> (\delta_i)^n is used but for fixed i, the same force F_i is used
> irrespective of the value of n.
> (See also, J. Chem. Phys. vol. 117, p2449 (2002), the 3rd line from the
> bottom of the 2nd column on the first page.)
The force on the non+central boxes are also stored separately in an
arrya fshift. The displacement part of the position is thus taken care
> 2. How does GROMACS calculate the virial for charged virtual sites using
> the Ewald summation (or SPME)? One needs to redistribute the forces
> before calculating the virial but if one uses this approach one
> will need not only 27 possible \delta_i but, technically, an infinite
The order of redistributing forces and virial calculation is important,
and one should indeed do it in the order you describe. This is not
completely obvious, since the position used for computing the forces are
those of the virtual sites. I do not see why one would need to compute
multiple deltas, but it could be that the contribution to fshift (see
above) is computed before redistributing the forces, which would not be
consistent in that case. GROMACS does reproduce e.g. the density of
TIP4P and TIP5P water correctly, so it seems that there is not a large
error if there is one at all.
> Thanks for your input.
> gmx-users mailing list gmx-users at gromacs.org
> 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
More information about the gromacs.org_gmx-users