I am trying to get the potential energy of a system at each frame of a molecular dynamics trajectory after running the simulation. My current approach is as follows:
- Use MDAnalysis to iterate through the trajectory frames and write out the atom coordinates to separate topology files
- Iterate through each of the files and create a simulation using OpenMM and then use their getPotentialEnergy() method to extract the energy from the simulation state
The issue is that when I test this using the same forcefield that was used for the simulation, I get potential energy values that differ on the order of 10% from the potential energy reported during the simulation. I have also noticed that the first topology file created by step one has atom coordinates that differ significantly from the input topology file, but I do not know why.
What could the issues with my implementation be, or is there an alternate way to accomplish this accurately?