symplectic integrator of any order

[Yucheng-Zhang]

1. We adapt the code for the flexibility of using a symplectic integrator of any order. The user could set the coefficients (e.g. the default [[0, 0, 1], [0, 0.5, 1]] for leapfrog) in configuration.
2. The loop over a_nbody is also switched from (a_prev, a_next) to index i. This is not necessary for the symplectic upgrade, but would make things easier for implementing observables.

[eelregit]

Adding to Yucheng's 2nd point, using i instead of a's in the time integration makes it easy to interpolate a snapshot between 2 time steps. No extra memory is needed to save both states for the interpolation. This trick won't work for light cones though, because the interception time needs to be solved using both states before the interpolation.

[eelregit]

Very useful reference: https://www.unige.ch/~hairer/poly_geoint/week2.pdf

[eelregit]

Actually, now I think snapshot interpolation is possible even with the observe(a_prev, a_next, ...) interface, which seems somewhat more modular to me.