A Geodesic Integrator

We may consider alternative splittings of the dynamics in order to build a library of available schemes, just as we did in the unconstrained case in Sect. 7.3.1. Solving the constrained OU part exactly would allieviate the condition (7.46), while we may hope to reduce discretization error in the spirit of the success of the BAOAB scheme in Sect. 7.9.3. The natural splitting template to consider would be to again split into potential and kinetic energies, with the constrained Ornstein-Uhlenbeck piece  

Observe that the holonomic constraint (as well as the hidden constraint) is maintained in each part, changing the usual A piece into a geodesic drift along the 

We then ask how to analyze such schemes using the algebraic machinery developed throughout this chapter. As we solve each vector field on the constraint manifold, it is natural to consider the associated Lie derivatives on these manifolds and solve the corresponding Fokker-Planck equation for the discretization to find an invariant distribution, just as we did in the unconstrained case. This is a complicated programme and we do not develop this in detail here. 

It is clear that the c- BAOAB] scheme gives a superior calculation of the average of potential energy (the computed distributions are plotted in Fig. 7.12), with an apparent order of magnitude improvement in the error compared to c-[OBABO. Additionally it allows for a significant improvement in stability between the two 

---