  We present heyoka, a new, modern and general-purpose implementation of
Taylor's integration method for the numerical solution of ordinary differential
equations. Detailed numerical tests focused on difficult high-precision
gravitational problems in astrodynamics and celestial mechanics show how our
general-purpose integrator is competitive with and often superior to
state-of-the-art specialised symplectic and non-symplectic integrators in both
speed and accuracy. In particular, we show how Taylor methods are capable of
satisfying Brouwer's law for the conservation of energy in long-term
integrations of planetary systems over billions of dynamical timescales. We
also show how close encounters are modelled accurately during simulations of
the formation of the Kirkwood gaps and of Apophis' 2029 close encounter with
the Earth (where heyoka surpasses the speed and accuracy of domain-specific
methods). heyoka can be used from both C++ and Python, and it is publicly
available as an open-source project.
