Feature #2288
gmx msd doesn't optimally handle missing input trajectory frames
Description
I am reporting an issue with the gmx msd command from gromacs version 5.1.2. The issue is that when the input .xtc file is missing a frame, the output skips that timestep (delta time) for the MSD. Such skipping of a delta time value corresponding to a missing time value is not necessary. I am not sure if it is an output problem or an assumption that the input time-series is continuous, or perhaps if everything is computed correctly but the output simply skips a timestep unnecessarily.
As I show below, this particular input .xtc file is missing a frame at time 261000:
Command line: gmx check -f TRAJ/bot201020ps_bot216020ps_vbot201020ps.xtc Checking file TRAJ/bot201020ps_bot216020ps_vbot201020ps.xtc Reading frame 0 time 0.000 # Atoms 197209 Precision 0.001 (nm) Reading frame 200 time 200000.000 Timesteps at t=260000 don't match (1000, 2000) Timesteps at t=262000 don't match (2000, 1000) Reading frame 1000 time 1001000.000
And here's evidence f the problem with the output: gmx msd with the command below thinks time follows the same numerical order as delta_time:
gmx msd -f $xtc -s $tpr -o out.xvg -ngroup 2 -lateral z -n index.ndx $ cat out.xvg ... 256000 6.14776 11.9849 257000 6.15575 12.0059 258000 6.1544 12.0186 259000 6.15638 12.0245 260000 6.1549 12.031 262000 6.15645 12.0408 <---- jump from 260000 to 262000 263000 6.15437 12.0461 264000 6.15383 12.0579 265000 6.14296 12.0726 ...
If the response it "garbage in, garbage out", then I suggest gmx msd throws an error and refuses to produce output unless things really are computed correctly.
History
#1 Updated by Berk Hess about 2 years ago
- Status changed from New to Accepted
the gmx_msd code is very difficult to understand, mainly because there is almost no documentation. But I think it computes MSD based on frame index offsets, so the results will be incorrect with non-uniform frame spacing.