Systematic slight discrepancy of the temperature value in NPT simulation

GROMACS version: 2020.3
GROMACS modification: No

Using Gromacs 2020.3 I noticed this anomaly in an NPT simulation of an aqueous solution (OPLS-AA force field; box size about 9x9x9 nm^3); mdrun runs on a Linux cluster using CPU + GPU). The average value of the temperature (controlled by the v-rescale thermostat) is slightly different from the value of the mdp directive (for example 302.77 K instead of 303 K, respectively). This systematic error is reproducible.
I found this behavior rather strange, since in past simulations the same thermostat setting gave the desired temperature value (in K) within the second or even the third decimal place.

These are the mdp directives I used:

; Temperature coupling
Tcoupl = v-rescale
nsttcouple = -1
nh-chain-length = 10
print-nose-hoover-chain-variables = no
; Groups to couple separately
tc_grps = System
; Time constant (ps) and reference temperature (K)
tau_t = 0.1
ref_t = 303
; pressure coupling
pcoupl = Parrinello-Rahman
pcoupltype = isotropic
nstpcouple = -1
; Time constant (ps), compressibility (1/bar) and reference P (bar)
tau_p = 2
compressibility = 4.4e-5
ref-p = 0.987
; Scaling of reference coordinates, No, All or COM
refcoord-scaling = all

and these the output from gmx energy on four subsequent trajectory intervals each of 100 ns

Energy Average Err.Est. RMSD Tot-Drift

Temperature 302.77 0.0014 1.14204 0.000128223 (K)

Temperature 302.768 0.0039 1.14271 0.00180873 (K)

Temperature 302.766 0.0025 1.14083 -0.011248 (K)

Temperature 302.764 0.0013 1.14201 -0.00211153 (K)

Obviously for my purposes to have 303 or 302.77 K does not make difference, however I would like to be sure that this is not the clue of possible problems in the simulation setup.

Did anyone detect a similar behaviour?

Thanks in advance for any comment or suggestion


Did you use exactly the same mdp options in both cases?
Which integrator did you use?
Could you compare the drift in the conserved energy quantity?