Polyhedral Potential and Variational Integrator Computation of the Full   Two Body Problem
  We present a combination of tools which allows for investigation of the coupled orbital and rotational dynamics of two rigid bodies with nearly arbitrary shape and mass distribution, under the influence of their mutual gravitational potential. Methods for calculating that mutual potential and resulting forces and moments for a polyhedral body representation are simple and efficient. Discrete equations of motion, referred to as the Lie Group Variational Integrator (LGVI), preserve the structure of the configuration space, SE(3), as well as the geometric features represented by the total energy and the total angular momentum. The synthesis of these approaches allows us to simulate the full two body problem accurately and efficiently. Simulation results are given for two octahedral rigid bodies for comparison with other integration methods and to show the qualities of the results thus obtained. A significant improvement is seen over other integration methods while correctly capturing the interesting effects of strong orbit and attitude dynamics coupling, in multiple scenarios. 
Eugene G. Fahnestock, Taeyoung Lee, Melvin Leok, N. Harris McClamroch,   Daniel J. Scheeres