Variational collision integrator for polymer chains

The numerical simulation of many-particle systems (e.g. in molecular dynamics) often involves constraints of various forms. We present a symplectic integrator for mechanical systems with holonomic (bilateral) and unilateral contact constraints, the latter being in the form of a non-penetration condition. The scheme is based on a discrete variant of Hamilton's principle in which both the discrete trajectory and the unknown collision time are varied (cf. [R. Fetecau, J. Marsden, M. Ortiz, M. West, Nonsmooth Lagrangian mechanics and variational collision integrators, SIAM J. Appl. Dyn. Syst. 2 (2003) 381-416]). As a consequence, the collision event enters the discrete equations of motion as an unknown that has to be computed on-the-fly whenever a collision is imminent. The additional bilateral constraints are efficiently dealt with employing a discrete null space reduction (including a projection and a local reparametrisation step) which considerably reduces the number of unknowns and improves the condition number during each time-step as compared to a standard treatment with Lagrange multipliers. We illustrate the numerical scheme with a simple example from polymer dynamics, a linear chain of beads, and test it against other standard numerical schemes for collision problems.