[gmx-users] Langevin sd integrator: wrong average temperature

[Andrey Frolov]

Dear GMX-users,

I perform MD simulation of one point on the phase space of supercritical CO2. I use NVT ensemble with the reference T=308.15K, number molecules 548, the box size is 3.6^3 nm^3. I use different thermostats: Berendsen, v-rescale, Nose-Hoover and Langevin (sd). The problem is that "sd" integrator gives a slightly elevated temperature (by ~ 3K) compared to the reference one. Also, the potential energy is higher for simulations with "sd" integrator. This strongly depends on the integration timestep, the smaller the timestep the lower the deviation. 

I would appreciate very much if someone gives a hint why "sd" behaves unexpectedly in this case ?
I should say that for the simulations of bulk tip3p water at ambient conditions I observe similar behavior (data not shown): the temperature tends to be slightly elevated with "sd" integrator (less than by ~ 1K, comparable to the statistical error)
               
Total simulation time for each parameter set is 5 ns, the averages for the last 4 ns:
                        
integrator,tcoupl,tau_t[ps],dt[ps]        Average   Err.Est.  RMSD  Tot-Drift
-------------------------------------------------------------------------------
md,ber,2.0,dt=0.002                          Temperature    308.2      0.0      3.4      0.0      (K)
md,nose,2.0,dt=0.002                       Temperature    308.1      0.0      6.1      0.0      (K)
md,vrescale,2.0,dt=0.002                 Temperature    308.1      0.1      5.9      0.6      (K)
sd,no,2,dt=0.0005                             Temperature    308.4      0.2      6.2     -0.3      (K)
sd,no,2,dt=0.001                               Temperature    308.8      0.1      6.2     -0.4      (K)
sd,no,2,dt=0.002                               Temperature    311.5      0.1      6.1      0.6      (K)


gromacs_4.6.5 (double precision). 
Also, I should note that the compressibility of the fluid is about an order higher than that of water. (~ 7.8e-4 bar^-1).

Please, see the example "*.tpr" files in a separate e-mail.

Please, also note that this is question is a continuation of my previous post and bug report on similar issue (which was not resolved), where I wrongly guessed that the problem might be related to extended ensemble barostat's degrees of freedom not substracted from the calculation of temperature. In this case no barostat applied, but the problem still remains.
http://comments.gmane.org/gmane.science.biology.gromacs.user/65935
http://redmine.gromacs.org/issues/1418


Thank you very much.
Andrey Frolov
----
Junior researcher 
Institute of solution chemistry 
Russian Academy of Sciences