Closed ChiCheng45 closed 1 month ago
Actually, I had a think and I realized that if you used something like the velocity verlet algorithm. The velocity should be equal to a central difference of the positions. This is probably why order=2 is the same as the velocity data from the MD simulation.
Differences are probably expected. However, it is interesting to note that perhaps in some situations it might be better to recalculate the velocities with a higher order method even when MD velocity data is available assuming all frames from the MD simulation are used.
Description of the error I've been testing the accuracies of the numerical derivatives with the cp2k water trajectory, this trajectory has velocity data so we can compare results. The timestep is 0.5 fs. Here is what I get for the temperatures.
temperature_0: order=0 temperature_2: order=2 temperature_3: order=3 temperature_4: order=4 temperature_5: order=5
order 1 and order 2 are the same since #369.
The odd thing is that order > 2 always gives slightly larger temperatures, could this be a bug in the code? The most accurate is order=1 or 2. Note that the first (and also the last) frame has to be calculated with forward numerical derivatives, this gives larger errors in the temperatures with order=1 or 2.
Suggested fix Replace all numerical derivatives so that it uses numpy.gradient which is essentially order=2 in MDANSE doing this will also resolve #369 since we fix it to order=2. Maybe also drop the first and last frames, although this is probably unnecessary since the error will disappear when we time average something.