Energy conservation for water

The constraint force algorithm has been used to integrate the equations of motion of 216 SPC/E[17] water molecules at ambient conditions with a time step $\delta t = 2$ fs. Results for a 20 ps interval are displayed in Fig. 1 as the ragged line with circles. The quantity $\langle E \rangle$ is the average energy for the 20 ps interval and $\Delta E = E(t)-\langle E \rangle$. The sloping line with squares is for the same system using a "scaling algorithm" where a velocity Verlet algorithm is used to integrate $Q_\alpha$ and $\dot Q_{\alpha}$ with scaling to impose the constraints. Clearly this demonstrates the superiority of the constraint force algorithm over the scaling algorithm. Note that a similar figure is given by Omelyan [10]. However in this case the running average $\langle E(t)\rangle$ instead of $\langle E \rangle$ is used in the denominator. Use of the running average can be misleading because it can mask a systematic drift in energy.

FIG. 1. The departures from the average energy for SPC/E water over a 20-ps time interval are shown for the algorithm discussed above as a line with circles, and for an algorithm where the constaints on the quaternions are imposed by scaling as a line with squares.