GROMACS

Is it center of mass of the whole system? That is, is the box changing shape?

Have you tried 4.6.3?

Yes, the center of mass of the whole system (which is a polymer chain with 400 monomer) is drifting. The box shape stays cubic and with a bit oscillation due to pressure control (The density is constant with correct value for polymer chain at all specified temperatures.).

I did simulations with 4.6 on a cluster and the only version available was 4.6 and below. If you think that may cause problem I can request an upgrade and repeat the calculations.

What's the shortest time scale over which it is statistically significant?

ABolfazl Noorjahan wrote:

Yes, the center of mass of the whole system (which is a polymer chain with 400 monomer) is drifting. The box shape stays cubic and with a bit oscillation due to pressure control (The density is constant with correct value for polymer chain at all specified temperatures.).

I did simulations with 4.6 on a cluster and the only version available was 4.6 and below. If you think that may cause problem I can request an upgrade and repeat the calculations.

Is your molecule cyclic?

David van der Spoel wrote:

Is your molecule cyclic?

It's a single polyethylene chain. I'm not sure what to make of that, other than perhaps something odd might happen with temperature groups. Note that it's using h-bonds, not all bonds.

Currently, it looks like an abrupt 1 nm vdw cutoff, which could give artifacts.

I'd like to see it run with (and I am trying a short version of this now).

nstcalctcoupl = 1
nstcalpcoupl = 1
nstcalcenergy = 1
nstcomm = 1

To remove any weird issues with properties being calculated only every N steps.

I started running this but realized that it may take a rather long time to check this -- it's running on my laptop at 2 ns/day, and it's not clear the issues would really be distinguishable for 1-2 days. Error cases that can be distinguished are usually better. . .

So far (20-30 ps in), it appears that the COM velocities are about 10x times larger for md-vv, using nstcomm = 100. For md/pr, it's down at 5^10-6, while for md-vv, it's more like 5^10-5, occasionally 1^10-4. This seems somewhat suspicious. I'll keep investigating.

Interestingly, it seems to be happening when nstcomm = 1 as well, and also when MTTK is off. If it's just that it's outputting the velocities at a different place in the loop, then that's not a big deal. I'll let it go overnight and see if the drift continues, and investigate the code further.