Symplectic Numerical Integrators in Constrained Hamiltonian Systems

Abstract

Recent work reported in the literature suggests that for the long-time integration of Hamiltonian dynamical systems one should use methods that preserve the symplectic (or canonical) structure of the flow. Here we investigate the symplecticness of numerical integrators for constrained dynamics, such as occur in molecular dynamics when bond lengths are made rigid in order to overcome stepsize limitations due to the highest frequencies. This leads to a constrained Hamiltonian system of smaller dimension. Previous work has shown that it is possible to have methods which are symplectic on the constraint manifold in phase space. Here it is shown that the very popular Verlet method with SHAKE-type constraints is equivalent to the same method with RATTLE-type constraints and that the latter is symplectic and time reversible. (This assumes that the iteration is carried to convergence.) We also demonstrate the global convergence of the Verlet scheme in the presence of SHAKE-type and RATTLE-type constraints. We perform numerical experiments to compare these methods with the second-order backward differentiation method, commonly recommended for ordinary differential equations with constraints.