[lammps-users] [Q] NVT integrator

Dear Steve,

I think I found a bug in velocity Verlet integrator in fix_nvt.cpp file.
In that file, the velocity is updated as

void FixNVT::final_integrate()
  for (int i = 0; i < nlocal; i++) {
    if (mask[i] & groupbit) {
      dtfm = dtf / mass[type[i]] * factor;
      v[i][0] = v[i][0]*factor + dtfm*f[i][0];
      v[i][1] = v[i][1]*factor + dtfm*f[i][1];
      v[i][2] = v[i][2]*factor + dtfm*f[i][2];

Instead, it should be updated as

      v[i][0] = (v[i][0] + dtfm*f[i][0])*factor;
      v[i][1] = (v[i][1] + dtfm*f[i][1])*factor;
      v[i][2] = (v[i][2] + dtfm*f[i][2])*factor;

I checked the reversibility and also found the same pattern in
fix_nph.cpp and fix_npt.cpp files.

Happy New Year,

Sorry Steve,

It's my bad.
The change that I've mentioned doesn't improve the reversibility of
Verlet integrator of NVT case at all.
Still I wonder why the reversibility of integrator for NVT isn't as
good as that for NVE.