  Multibody dynamics simulators are an important tool in many fields, including
learning and control for robotics. However, many existing dynamics simulators
suffer from inaccuracies when dealing with constrained mechanical systems due
to unsuitable integrators with bad energy behavior and problematic constraint
violations, for example for contact interactions. Variational integrators are
numerical discretization methods that can reduce physical inaccuracies when
simulating mechanical systems, and formulating the dynamics in maximal
coordinates allows for easy and numerically robust incorporation of constraints
such as kinematic loops or contacts. Therefore, this article derives a
variational integrator for mechanical systems with equality and inequality
constraints in maximal coordinates. Additionally, efficient graph-based
sparsity-exploiting algorithms for solving the integrator are provided and
implemented as an open-source simulator. The evaluation of the simulator shows
improved physical accuracy due to the variational integrator and the advantages
of the sparse solvers. Comparisons to minimal-coordinate algorithms show
improved numerical robustness and application examples of a walking robot and
an exoskeleton with explicit constraints demonstrate the necessity and
capabilities of maximal coordinates.
