- Generated by
1.9.4
|
GeographicLib 2.1.2
|
The NormalGravity class computes "normal gravity" which refers to the exact (classical) solution of an idealised system consisting of an ellipsoid of revolution with the following properties:
(N.B. The mass always appears in the combination GM, with units m3/s2, where G is the gravtitational constant.) The distribution of the mass M within the ellipsoid is such that the surface of the ellipsoid is at a constant normal potential where the normal potential is the sum of the gravitational potential (due to the gravitational attraction) and the centrifugal potential (due to the rotation). The resulting field exterior to the ellipsoid is called normal gravity and was found by Somigliana (1929). Because the potential is constant on the ellipsoid it is sometimes referred to as the level ellipsoid.
Go to
References:
Two set of formulas are presented: those of Heiskanen and Moritz (1967) which are applicable to an oblate ellipsoid and a second set where the variables are distinguished with primes which apply to a prolate ellipsoid. The primes are omitted from those variables which are the same in the two cases. In the text, the parenthetical "alt." clauses apply to prolate ellipsoids.
Cylindrical coordinates \( R,Z \) are expressed in terms of ellipsoidal coordinates
\[ \begin{align} R &= \sqrt{u^2 + E^2} \cos\beta = u' \cos\beta,\\ Z &= u \sin\beta = \sqrt{u'^2 + E'^2} \sin\beta, \end{align} \]
where
\[ \begin{align} E^2 = a^2 - b^2,&{} \qquad u^2 = u'^2 + E'^2,\\ E'^2 = b^2 - a^2,&{} \qquad u'^2 = u^2 + E^2. \end{align} \]
Surfaces of constant \( u \) (or \( u' \)) are confocal ellipsoids. The surface \( u = 0 \) (alt. \( u' = 0 \)) corresponds to the focal disc of diameter \( 2E \) (alt. focal rod of length \( 2E' \)). The level ellipsoid is given by \( u = b \) (alt. \( u' = a \)). On the level ellipsoid, \( \beta \) is the familiar parametric latitude.
In writing the potential and the gravity, it is useful to introduce the functions
\[ \begin{align} Q(z) &= \frac1{2z^3} \biggl[\biggl(1 + \frac3{z^2}\biggr)\tan^{-1}z - \frac3z\biggr],\\ Q'(z') &= \frac{(1+z'^2)^{3/2}}{2z'^3} \biggl[\biggl(2 + \frac3{z'^2}\biggr)\sinh^{-1}z' - \frac{3\sqrt{1+z'^2}}{z'}\biggr],\\ H(z) &= \biggl(3Q(z)+z\frac{dQ(z)}{dz}\biggr)(1+z^2)\\ &= \frac1{z^4}\biggl[3(1+z^2) \biggl(1-\frac{\tan^{-1}z}z\biggr)-z^2\biggr],\\ H'(z') &= \frac{3Q'(z')}{1+z'^2}+z'\frac{dQ'(z')}{dz'}\\ &= \frac{1+z'^2}{z'^4} \biggl[3\biggl(1-\frac{\sqrt{1+z'^2}\sinh^{-1}z'}{z'}\biggr) +z'^2\biggr]. \end{align} \]
The function arguments are \( z = E/u \) (alt. \( z' = E'/u' \)) and the unprimed and primed quantities are related by
\[ \begin{align} Q'(z') = Q(z),&{} \qquad H'(z') = H(z),\\ z'^2 = -\frac{z^2}{1 + z^2},&{} \qquad z^2 = -\frac{z'^2}{1 + z'^2}. \end{align} \]
The functions \( q(u) \) and \( q'(u) \) used by Heiskanen and Moritz are given by
\[ q(u) = \frac{E^3}{u^3}Q\biggl(\frac Eu\biggr),\qquad q'(u) = \frac{E^2}{u^2}H\biggl(\frac Eu\biggr). \]
The functions \( Q(z) \), \( Q'(z') \), \( H(z) \), and \( H'(z') \) are more convenient for use in numerical codes because they are finite in the spherical limit \( E \rightarrow 0 \), i.e., \( Q(0) = Q'(0) = \frac2{15} \), and \( H(0) = H'(0) = \frac25 \).
The normal potential is the sum of three components, a mass term, a quadrupole term and a centrifugal term,
\[ U = U_m + U_q + U_r. \]
The mass term is
\[ U_m = \frac {GM}E \tan^{-1}\frac Eu = \frac {GM}{E'} \sinh^{-1}\frac{E'}{u'}; \]
the quadrupole term is
\[ \begin{align} U_q &= \frac{\omega^2}2 \frac{a^2 b^3}{u^3} \frac{Q(E/u)}{Q(E/b)} \biggl(\sin^2\beta-\frac13\biggr)\\ &= \frac{\omega^2}2 \frac{a^2 b^3}{(u'^2+E'^2)^{3/2}} \frac{Q'(E'/u')}{Q'(E'/a)} \biggl(\sin^2\beta-\frac13\biggr); \end{align} \]
finally, the rotational term is
\[ U_r = \frac{\omega^2}2 R^2 = \frac{\omega^2}2 (u^2 + E^2) \cos^2\beta = \frac{\omega^2}2 u'^2 \cos^2\beta. \]
\( U_m \) and \( U_q \) are both gravitational potentials (due to mass within the ellipsoid). The total mass contributing to \( U_m \) is \( M \); the total mass contributing to \( U_q \) vanishes (far from the ellipsoid, the \( U_q \) decays inversely as the cube of the distance).
\( U_m \) and \( U_q + U_r \) are separately both constant on the level ellipsoid. \( U_m \) is the normal potential for a massive non-rotating ellipsoid. \( U_q + U_r \) is the potential for a massless rotating ellipsoid. Combining all the terms, \( U \) is the normal potential for a massive rotating ellipsoid.
The figure shows normal gravity for the case \( GM = 1 \), \( a = 1 \), \( b = 0.8 \), \( \omega = 0.3 \). The level ellipsoid is shown in red. Contours of constant gravity potential are shown in blue; the contour spacing is constant outside the ellipsoid and equal to 1/20 of the difference between the potentials on the ellipsoid and at the geostationary point ( \( R = 2.2536 \), \( Z = 0 \)); inside the ellipsoid the contour spacing is 5 times greater. The green lines are stream lines for the gravity; these are spaced at intervals of 10° in parametric latitude on the surface of the ellipsoid. The normal gravity is continued into the level ellipsoid under the assumption that the mass is concentrated on the focal disc, shown in black.
Typically, the normal potential, \( U \), is only of interest for outside the ellipsoid \( u \ge b \) (alt. \( u' \ge a \)). In planetary applications, an open problem is finding a mass distribution which is in hydrostatic equilibrium (the mass density is non-negative and a non-decreasing function of the potential interior to the ellipsoid).
However it is possible to give singular mass distributions consistent with the normal potential.
For a non-rotating body, the potential \( U = U_m \) is generated by a sandwiching the mass \( M \) uniformly between the level ellipsoid with semi-axes \( a \) and \( b \) and a close similar ellipsoid with semi-axes \( (1-\epsilon)a \) and \( (1-\epsilon)b \). Chasles (1840) extends a theorem of Newton to show that the field interior to such an ellipsoidal shell vanishes. Thus the potential on the ellipsoid is constant, i.e., it is indeed a level ellipsoid. This result also holds for a non-rotating triaxial ellipsoid.
Observing that \( U_m \) and \( U_q \) as defined above obey \( \nabla^2 U_m = \nabla^2 U_q = 0 \) everywhere for \( u > 0 \) (alt. \( u' > 0 \)), we see that these potentials correspond to masses concentrated at \( u = 0 \) (alt. \( u' = 0 \)).
In the oblate case, \( U_m \) is generated by a massive disc at \( Z = 0 \), \( R < E \), with mass density (mass per unit area) \( \rho_m \) and moments of inertia about the equatorial (resp. polar) axis of \( A_m \) (resp. \( C_m \)) given by
\[ \begin{align} \rho_m &= \frac M{2\pi E\sqrt{E^2 - R^2}},\\ A_m &= \frac {ME^2}3, \\ C_m &= \frac {2ME^2}3, \\ C_m-A_m &= \frac {ME^2}3. \end{align} \]
This mass distribution is the same as that produced by projecting a uniform spherical shell of mass \( M \) and radius \( E \) onto the equatorial plane.
In the prolate case, \( U_m \) is generated by a massive rod at \( R = 0 \), \( Z < E' \) and now the mass density \( \rho'_m \) has units mass per unit length,
\[ \begin{align} \rho'_m &= \frac M{2E'},\\ A_m &= \frac {ME'^2}3, \\ C_m &= 0, \\ C_m-A_m &= -\frac {ME'^2}3. \end{align} \]
This mass distribution is the same as that produced by projecting a uniform spherical shell of mass \( M \) and radius \( E' \) onto the polar axis.
Similarly, \( U_q \) is generated in the oblate case by
\[ \begin{align} \rho_q &= \frac{a^2 b^3 \omega^2}G \frac{2E^2 - 3R^2}{6\pi E^5 \sqrt{E^2 - R^2} Q(E/b)}, \\ A_q &= -\frac{a^2 b^3 \omega^2}G \frac2{45\,Q(E/b)}, \\ C_q &= -\frac{a^2 b^3 \omega^2}G \frac4{45\,Q(E/b)}, \\ C_q-A_q &= -\frac{a^2 b^3 \omega^2}G \frac2{45\,Q(E/b)}. \end{align} \]
The corresponding results for a prolate ellipsoid are
\[ \begin{align} \rho_q' &= \frac{a^2 b^3 \omega^2}G \frac{3Z^2 - E'^2}{12\,E'^5 Q'(E'/a)}, \\ A_q &= \frac{a^2 b^3 \omega^2}G \frac2{45\,Q'(E'/a)}, \\ C_q &= 0, \\ C_q-A_q &= -\frac{a^2 b^3 \omega^2}G \frac2{45\,Q'(E'/a)}. \end{align} \]
Summing up the mass and quadrupole terms, we have
\[ \begin{align} A &= A_m + A_q, \\ C &= C_m + C_q, \\ J_2 & = \frac{C - A}{Ma^2}, \end{align} \]
where \( J_2 \) is the dynamical form factor.
Each term in the potential contributes to the gravity on the surface of the ellipsoid
\[ \gamma = \gamma_m + \gamma_q + \gamma_r; \]
These are the components of gravity normal to the ellipsoid and, by convention, \( \gamma \) is positive downwards. The tangential components of the total gravity and that due to \( U_m \) vanish. Those tangential components of the gravity due to \( U_q \) and \( U_r \) cancel one another.
The gravity \( \gamma \) has the following dependence on latitude
\[ \begin{align} \gamma &= \frac{b\gamma_a\cos^2\beta + a\gamma_b\sin^2\beta} {\sqrt{b^2\cos^2\beta + a^2\sin^2\beta}}\\ &= \frac{a\gamma_a\cos^2\phi + b\gamma_b\sin^2\phi} {\sqrt{a^2\cos^2\phi + b^2\sin^2\phi}}, \end{align} \]
and the individual components, \( \gamma_m \), \( \gamma_q \), and \( \gamma_r \), have the same dependence on latitude. The equatorial and polar gravities are
\[ \begin{align} \gamma_a &= \gamma_{ma} + \gamma_{qa} + \gamma_{ra},\\ \gamma_b &= \gamma_{mb} + \gamma_{qb} + \gamma_{rb}, \end{align} \]
where
\[ \begin{align} \gamma_{ma} &= \frac{GM}{ab},\qquad \gamma_{mb} = \frac{GM}{a^2},\\ \gamma_{qa} &= -\frac{\omega^2 a}6 \frac{H(E/b)}{Q(E/b)} = -\frac{\omega^2 a}6 \frac{H'(E'/a)}{Q'(E'/a)},\\ \gamma_{qb} &= \frac{\omega^2 b}3 \frac{H(E/b)}{Q(E/b)} = \frac{\omega^2 b}3 \frac{H'(E'/a)}{Q'(E'/a)},\\ \gamma_{ra} &= -\omega^2 a,\qquad \gamma_{rb} = 0. \end{align} \]
Performing an average of the surface gravity over the area of the ellipsoid gives
\[ \langle \gamma \rangle = \frac {4\pi a^2 b}A \biggl(\frac{2\gamma_a}{3a} + \frac{\gamma_b}{3b}\biggr), \]
where \( A \) is the area of the ellipsoid
\[ \begin{align} A &= 2\pi\biggl( a^2 + ab\frac{\sinh^{-1}(E/b)}{E/b} \biggr)\\ &= 2\pi\biggl( a^2 + b^2\frac{\tan^{-1}(E'/a)}{E'/a} \biggr). \end{align} \]
The contributions to the mean gravity are
\[ \begin{align} \langle \gamma_m \rangle &= \frac{4\pi}A GM, \\ \langle \gamma_q \rangle &= 0 \quad \text{(as expected)}, \\ \langle \gamma_r \rangle &= -\frac{4\pi}A \frac{2\omega^2 a^2b}3,\\ \end{align} \]
resulting in
\[ \langle \gamma \rangle = \frac{4\pi}A \biggl(GM - \frac{2\omega^2 a^2b}3\biggr). \]
The solution for the normal gravity is well defined for arbitrary \( M \), \( \omega \), \( a > 0\), and \( f < 1 \). (Note that arbitrary oblate and prolate ellipsoids are possible, although hydrostatic equilibrium would not result in a prolate ellipsoid.) However, it is much easier to measure the dynamical form factor \( J_2 \) (from the motion of artificial satellites) than the flattening \( f \). (Note too that \( GM \) is typically measured from satellite or astronomical observations and so it includes the mass of the atmosphere.)
So a question for the software developer is: given values of \( M > 0\), \( \omega \), and \( a > 0 \), what are the allowed values of \( J_2 \)? We restrict the question to \( M > 0 \). The (unphysical) case \( M = 0 \) is problematic because \( M \) appears in the denominator in the definition of \( J_2 \). In the (also unphysical) case \( M < 0 \), a given \( J_2 \) can result from two distinct values of \( f \).
Holding \( M > 0\), \( \omega \), and \( a > 0 \) fixed and varying \( f \) from \( -\infty \) to \( 1 \), we find that \( J_2 \) monotonically increases from \( -\infty \) to
\[ \frac13 - \frac8{45\pi} \frac{\omega^2 a^3}{GM}. \]
Thus any value of \( J_2 \) less that this value is permissible (but some of these values may be unphysical). In obtaining this limiting value, we used the result \( Q(z \rightarrow \infty) \rightarrow \pi/(4 z^3) \). The value
\[ J_2 = -\frac13 \frac{\omega^2 a^3}{GM} \]
results in a sphere ( \( f = 0 \)).
1.9.4