Add normal_gravitational_potential, normal_gravity_potential, and cen…

…trifugal_potential (#187)

This commit makes several additions in order to compute the normal
gravity/gravitational potential on or above the ellipsoid.

* Adds `normal_gravity_potential()`, `normal_gravitational_potential()`
and `centrifugal_potential()` methods to both the Ellipsoid and Sphere
class.
* Adds two conversion routines, to and from geodetic and
ellipsoidal-harmonic coordinates:
`Ellipsoid.geodetic_to_ellipsoidal_harmonic()` and
`Ellipsoid.ellipsoidal_harmonic_to_geodetic()`.
* Adds tests to ensure that `normal_gravity_potential() =
normal_gravitational_potential() + centrifugal_potential()`.
* Adds test for the coordinate transforms. See below: these are
surprisingly not as accurate as I thought they would be.

The centrifugal potential calculation needs the perpendicular
distance to the rotation axis. For this, I used the definition of the
prime radius of curvature and geodetic latitude which gives: $(N(\phi) +
h) \cos\phi$, where $\phi$ is geodetic latitude and $h$ is ellipsoidal
height. There are probably other ways that this could be calculated.

The geodetic to ellipsoidal harmonic coordinate transforms work, but
I am worried about the precision. The conversion from geodetic to
ellipsoidal uses the same code as in the `normal_gravity` method, which
is based on the equations in Lakshmanan (1991). The inverse is just a a
couple atans for the latitude, and the ellipsoidal height is computed as
the difference of the `prime_vertical_radius` of the two ellipsoids.
Even if I can understand why the height might lose precision this way,
the latitude shouldn't.

Maybe this is good enough. I suspect that the loss of precision comes
from the Lakshmanan equations (which have some subtractions of similarly
sized numbers). Perhaps there is nothing we can do about it.
Nevertheless, I note that if this is a problem, then it will also affect
the normal gravity routine that uses the same method.