Energies are not obtained at the same timestep.

[kosukesato]

https://github.com/Quantum-Dynamics-Hub/Libra-X/blob/master/src/md.py I have a question about the part for obtaining energies (see 360-384 lines in the attached link). In the concerned code, kinetic energy (epot[cnt]) and potential energy (ekin[cnt]) are obtained at "t+dt/2", so the sum of them (etot[cnt]) is at "t+dt/2". But, thermostat energy(therm[cnt].energy()) is at "t+dt". We don't worry about that when we check energy conservation under NVE-MD. But, when we check under NVT-MD, we will care about the time difference between etot[cnt] and therm[cnt].energy(). To avoid this problem, the energies should be obtained at the same time. If you have any comments on this issue, please let me know. Thank you in advance.

P.S. Actually, I tried to open this issue in Chemistry Stack Exchange for testing the system, but I failed to do so in my company (maybe connection problem occured).

[alexvakimov]

Hi Kosuke, The potential energy is evaluated at q(t+dt), that is at the end of the time step, even though it is done in the middle of the integrator. I agree that the kinetic energy contribution comes at the time step t+dt/2. This should be easy to correct. Only, keep in mind (if you decide doing that), that we still need the kinetic energy at t+dt/2, because it is used to update thermostat "forces".

[kosukesato]

We use the potential energy at t+dt/2 because it is given by vibronic hamiltonian having the averaged one epot(t+dt/2) = 1/2 *(E(t)+E(t+dt)) via "compute_forces" function. I'll modify the part following your advice. Thank you.

[alexvakimov]

Thank you, Kosuke, for the modification!