fix timesteps to correct for the wrong size of the arrays

of previous_stardatas. The number stored is set by the macro
Solver_n_previous_stardatas(stardata) which can then depend on
e.g. preferences.

Added an RK2 solver : this should be considered experimental.
This works by calling evolution() twice with half timesteps.
The derivatives are applied using a new function apply_derivative
which uses some funky pointer magic to do the calculations.

Differences are tiny (3rd or 4th decimal places), while runtime
differences are currently large (x10) even in systems with lots of tides.

Events (SNe, comenvs) are not well handled at present, hence it is ONLY
experimental. It does work though :)