[gmx-users] Bugs in calculation of virial
Vincent Ballenegger
vcb25 at cam.ac.uk
Tue Jul 20 11:40:01 CEST 2004
Dear all,
I think I found three bugs in the way the virial (and hence the
pressure) is calculated in gromacs.
1) Bug in long range correction term:
The virial in gromacs is defined as (eq. 3.18 of manual)
tensor(Xi) = - (1/2) sum_{i<j} vect{r_ij} x vect{F_ij}.
The factor (1/2) has then be forgotten in equations (C.8) and (C.10) of
the manual, as well as in the code in sim_util.c that computes this LR
correction term (C.10).
2) Bug in long range correction term (again)
The correction (C.10) is to be applied to the scalar virial Xi. Gromacs
adds this term to each element on the diagonal of the virial tensor, so
a factor 1/3 is missing, since Xi = Trace(tensor(Xi)).
This bug was probably induced by the fact that the scalar pressure is
defined as *one third of* the trace of the pressure tensor, so it's
right to add the LR correction term (C.11) to each diagonal element of
the pressure tensor.
3) virial and external electric field
When a uniform external electric field is applied on the system, the
pressure should be computed from the *internal* virial, which excludes
the contribution from the forces due to the external electric field. As
far as I can see in the code, gromacs does not exclude this contribution
(which can be big if the applied field is strong).
In the revision notes, it is written that the virial contribution of
forces from an electric field has been removed in version 3.0.2. I
compiled however versions 3.0.1 and 3.0.2, and both give identical
results for the virial of SPC water in a strong electric field, so it
seems that this bug has not been corrected in version 3.0.2 for some
strange reason.
Regards,
Vincent Ballenegger
More information about the gromacs.org_gmx-users
mailing list