View source for Physics:Geopotential model

{{Short description|Theoretical description of Earth's gravimetric shape}}
In [[Physics:Geophysics|geophysics]] and [[Physics:Physical geodesy|physical geodesy]], a '''geopotential model''' is the theoretical analysis of measuring and calculating the effects of [[Astronomy:Earth|Earth]]'s [[Physics:Gravitational field|gravitational field]] (the [[Physics:Geopotential|geopotential]]).

==Newton's law==

[[Image:NewtonsLawOfUniversalGravitation.svg|thumb|Diagram of two masses attracting one another]]

[[Physics:Newton's law of universal gravitation|Newton's law of universal gravitation]] states that the gravitational force ''F'' acting between two point masses ''m''<sub>1</sub> and ''m''<sub>2</sub> with [[Physics:Centre of mass|centre of mass]] separation ''r'' is given by

:<math>\mathbf{F} = - G \frac{m_1 m_2}{r^2}\mathbf{\hat{r}}</math>

where ''G'' is the [[Physics:Gravitational constant|gravitational constant]] and '''r̂''' is the radial [[Unit vector|unit vector]]. For a non-pointlike object of continuous mass distribution, each mass element ''dm'' can be treated as mass distributed over a small volume, so the [[Volume integral|volume integral]] over the extent of object 2 gives:

{{NumBlk|:|<math>
\mathbf{\bar{F}} = - Gm_1 \int\limits_V \frac{\rho_2 }{r^2}\mathbf{\hat{r}}\,dx\,dy\,dz
</math>|{{EquationRef|1}}}}

with corresponding [[Physics:Gravitational potential|gravitational potential]]

{{NumBlk|:|<math>
u = - G \int\limits_V \frac{\rho_2}{r}\, dx\,dy\,dz
</math>|{{EquationRef|2}}}}

where ρ{{sub|2}} = ρ(''x'', ''y'', ''z'') is the [[Physics:Mass density|mass density]] at the [[Volume element|volume element]] and of the direction from the volume element to point mass 1. <math>u</math> is the gravitational potential energy per unit mass.

===The case of a homogeneous sphere===

In the special case of a sphere with a spherically symmetric mass density then ρ = ρ(''s''); i.e., density depends only on the radial distance

:<math>s = \sqrt{x^2 + y^2 + z^2} \,.</math>

These integrals can be evaluated analytically. This is the [[Physics:Shell theorem|shell theorem]] saying that in this case:

{{NumBlk|:|<math>
\bar{F} = -\frac{GMm}{R^2}\ \hat{r}
</math>|{{EquationRef|3}}}}

with corresponding [[Physics:Potential|potential]]

{{NumBlk|:|<math>
u = -\frac{GM}{r}
</math>|{{EquationRef|4}}}}

where ''M'' = ∫<sub>V</sub>ρ(''s'')''dxdydz'' is the total mass of the sphere.

==Spherical harmonics representation==

In reality, Earth is not exactly spherical, mainly because of its rotation around the polar axis that makes its shape slightly oblate. If this shape were perfectly known together with the exact mass density ρ = ρ(''x'', ''y'', ''z''), the integrals ({{EquationNote|1}}) and ({{EquationNote|2}}) could be evaluated with numerical methods to find a more accurate model for Earth's gravitational field. However, the situation is in fact the opposite. By observing the orbits of spacecraft and the Moon, Earth's gravitational field can be determined quite accurately and the best estimate of Earth's mass is obtained by dividing the product ''GM'' as determined from the analysis of spacecraft orbit with a value for ''G'' determined to a lower relative accuracy using other physical methods.

===Background===
{{main|Physics:Spherical harmonics}}
From the defining equations ({{EquationNote|1}}) and ({{EquationNote|2}}) it is clear (taking the partial derivatives of the integrand) that outside the body in empty space the following differential equations are valid for the field caused by the body:
{{NumBlk|:|<math>
\frac{\partial F_x }{\partial x} + \frac{\partial F_y }{\partial y} + \frac{\partial F_z }{\partial z} = 0
</math>|{{EquationRef|5}}}}
{{NumBlk|:|<math>
\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial z^2} = 0 </math>|{{EquationRef|6}}}}

Functions of the form <math>\phi = R(r)\, \Theta(\theta)\, \Phi(\varphi)</math> where (''r'', θ, φ) are the spherical coordinates which satisfy the partial differential equation ({{EquationNote|6}}) (the [[Laplace's equation|Laplace equation]]) are called spherical harmonic functions.

They take the forms:
{{NumBlk|:|<math>\begin{align}
  g(x, y, z) &= \frac{1}{r^{n+1}} P^m_n(\sin \theta) \cos m\varphi \,,&
    0 &\le m \le n \,,&
    n &= 0, 1, 2, \dots \\
  h(x, y, z) &= \frac{1}{r^{n+1}}  P^m_n(\sin \theta) \sin m\varphi \,,&
    1 &\le m \le n \,,&
    n &= 1, 2, \dots
\end{align}</math>|{{EquationRef|7}}}}

where spherical coordinates (''r'', θ, φ) are used, given here in terms of cartesian (''x, y, z'') for reference:
{{NumBlk|:|<math>\begin{align}
  x &= r \cos \theta \cos \varphi \\
  y &= r \cos \theta \sin \varphi \\
  z &= r \sin \theta\,,
\end{align}</math>|{{EquationRef|8}}}}

also ''P''<sup>0</sup><sub>''n''</sub> are the [[Legendre polynomials]] and ''P<sup>m</sup><sub>n</sub>'' for {{nowrap|1 ≤ ''m'' ≤ ''n''}} are the [[Physics:Associated Legendre polynomials|associated Legendre functions]].

The first spherical harmonics with ''n'' = 0, 1, 2, 3 are presented in the table below.

:{| class="wikitable"
|-
! ''n'' !! Spherical harmonics 
|-
| 0
| <math>\frac{1}{r}</math>
|-
|rowspan="3"| 1
| <math>\frac{1}{r^2} P^0_1(\sin\theta) = \frac{1}{r^2} \sin\theta</math>
|-
| <math>\frac{1}{r^2} P^1_1(\sin\theta) \cos\varphi= \frac{1}{r^2} \cos\theta \cos\varphi</math>
|-
| <math>\frac{1}{r^2} P^1_1(\sin\theta) \sin\varphi= \frac{1}{r^2} \cos\theta \sin\varphi</math>
|-
|rowspan="5"| 2
| <math>\frac{1}{r^3} P^0_2(\sin\theta) = \frac{1}{r^3} \frac{1}{2} (3\sin^2\theta - 1)</math>
|-
| <math>\frac{1}{r^3} P^1_2(\sin\theta) \cos\varphi = \frac{1}{r^3} 3 \sin\theta \cos\theta\ \cos\varphi</math>
|-
| <math>\frac{1}{r^3} P^1_2(\sin\theta) \sin\varphi = \frac{1}{r^3} 3 \sin\theta \cos\theta \sin\varphi</math>
|-
| <math>\frac{1}{r^3} P^2_2(\sin\theta) \cos2\varphi = \frac{1}{r^3} 3 \cos^2 \theta\ \cos2\varphi</math>
|-
| <math>\frac{1}{r^3} P^2_2(\sin\theta) \sin2\varphi = \frac{1}{r^3} 3 \cos^2 \theta \sin 2\varphi</math>
|-
|rowspan="7"| 3
| <math>\frac{1}{r^4} P^0_3(\sin\theta) = \frac{1}{r^4} \frac{1}{2} \sin\theta\ (5\sin^2\theta -3)</math>
|-
| <math>\frac{1}{r^4} P^1_3(\sin\theta) \cos\varphi = \frac{1}{r^4} \frac{3}{2}\ (5\ \sin^2\theta - 1) \cos\theta \cos\varphi</math>
|-
| <math>\frac{1}{r^4} P^1_3(\sin\theta) \sin\varphi = \frac{1}{r^4} \frac{3}{2}\ (5\ \sin^2\theta - 1) \cos\theta \sin\varphi</math>
|-
| <math>\frac{1}{r^4} P^2_3(\sin\theta) \cos 2\varphi = \frac{1}{r^4} 15 \sin\theta \cos^2 \theta \cos 2\varphi</math>
|-
| <math>\frac{1}{r^4} P^2_3(\sin\theta) \sin 2\varphi = \frac{1}{r^4} 15 \sin\theta \cos^2 \theta \sin 2\varphi</math>
|-
| <math>\frac{1}{r^4} P^3_3(\sin\theta) \cos 3\varphi = \frac{1}{r^4} 15 \cos^3 \theta \cos 3\varphi</math>
|-
| <math>\frac{1}{r^4} P^3_3(\sin\theta) \sin 3\varphi = \frac{1}{r^4} 15 \cos^3 \theta \sin 3\varphi</math>
|}

===Application===
The model for Earth's gravitational potential is a sum

{{NumBlk|:|<math>
u = -\frac{\mu}{r} + \sum_{n=2}^{N_z} \frac{J_n P^0_n(\sin\theta) }{r^{n+1}} + \sum_{n=2}^{N_t} \sum_{m=1}^n \frac{ P^m_n(\sin\theta) \left(C_n^m \cos m\varphi + S_n^m \sin m\varphi\right)}{r^{n+1}}
</math>|{{EquationRef|9}}}}

where <math>\mu = GM</math> and the coordinates ({{EquationNote|8}}) are relative to the standard geodetic reference system extended into space with origin in the center of the [[Earth:Reference ellipsoid|reference ellipsoid]] and with ''z''-axis in the direction of the polar axis.

The '''zonal terms''' refer to terms of the form:
:<math>\frac{P^0_n(\sin\theta)}{r^{n+1}} \quad n=0,1,2,\dots</math>

and the '''tesseral terms''' terms refer to terms of the form:
:<math>\frac{P^m_n(\sin\theta) \cos m\varphi}{r^{n+1}}\,, \quad 1 \le m \le n \quad n=1,2,\dots</math>
:<math>\frac{P^m_n(\sin\theta) \sin m\varphi}{r^{n+1}}</math>

The zonal and tesseral terms for ''n''&nbsp;=&nbsp;1 are left out in ({{EquationNote|9}}).  The coefficients for the n=1 with both m=0 and m=1 term correspond to an arbitrarily oriented dipole term in the multi-pole expansion. Gravity does not physically exhibit any dipole character and so the integral characterizing ''n''&nbsp;=&nbsp;1 must be zero.

The different coefficients ''J<sub>n</sub>'', ''C<sub>n</sub><sup>m</sup>'', ''S<sub>n</sub><sup>m</sup>'', are then given the values for which the best possible agreement between the computed and the observed spacecraft orbits is obtained.

As ''P''<sup>0</sup><sub>n</sub>(''x'') = −''P''<sup>0</sup><sub>n</sub>(−''x'') non-zero coefficients ''J<sub>n</sub>'' for '''odd ''n''''' correspond to a lack of symmetry "north–south" relative the equatorial plane for the mass distribution of Earth. Non-zero coefficients ''C<sub>n</sub><sup>m</sup>'', ''S<sub>n</sub><sup>m</sup>'' correspond to a lack of rotational symmetry around the polar axis for the mass distribution of Earth, i.e. to a "tri-axiality" of Earth.

For large values of ''n'' the coefficients above (that are divided by ''r''<sup>(''n'' + 1)</sup> in ({{EquationNote|9}})) take very large values when for example kilometers and seconds are used as units. In the literature it is common to introduce some arbitrary "reference radius" ''R'' close to Earth's radius and to work with the dimensionless coefficients

:<math>\begin{align}
      \tilde{J_n} &= -\frac{J_n}{\mu\ R^n}, &
  \tilde{C_{n}^m} &= -\frac{C_{n}^m}{\mu\ R^n}, &
  \tilde{S_{n}^m} &= -\frac{S_{n}^m}{\mu\ R^n}
\end{align}</math>

and to write the potential as

{{NumBlk|:|<math>
u = -\frac{\mu }{r} \left(1 + \sum_{n=2}^{N_z}  \frac{\tilde{J_n} P^0_n(\sin\theta) }{{(\frac{r}{R})}^n} + \sum_{n=2}^{N_t} \sum_{m=1}^n \frac{ P^m_n(\sin\theta) (\tilde{C_{n}^m} \cos m\varphi + \tilde{S_{n}^m} \sin m\varphi)}{{(\frac{r}{R})}^n}\right)
</math>|{{EquationRef|10}}}}

===Largest terms===
The dominating term (after the term −μ/''r'') in ({{EquationNote|9}}) is the "''J''<sub>2</sub> coefficient", representing the oblateness of Earth:

:<math>u = \frac{J_2\ P^0_2(\sin\theta)}{r^3} = J_2 \frac{1}{r^3} \frac{1}{2} (3\sin^2\theta -1) = J_2 \frac{1}{r^5} \frac{1}{2} (3 z^2 -r^2)</math>

Relative the coordinate system
{{NumBlk|:|<math>\begin{align}
  \hat{\varphi} &= -\sin \varphi \hat{x} + \cos \varphi \hat{y} \\
   \hat{\theta} &= -\sin \theta\, (\cos \varphi \hat{x} + \sin \varphi \hat{y}) + \cos\theta \hat{z} \\
        \hat{r} &=  \cos \theta\, (\cos \varphi \hat{x} + \sin \varphi \hat{y}) + \sin\theta \hat{z}
\end{align}</math>|{{EquationRef|11}}}}

[[File:Spherical coordinates unit vectors.svg|thumb|right|Figure 1: The unit vectors. This is wrong. There should be a theta, not lambda <math>\hat{\varphi}\ ,\ \hat{\theta}\ ,\ \hat{r}</math>]]

illustrated in figure 1 the components of the force caused by the "''J''<sub>2</sub> term" are
{{NumBlk|:|<math>\begin{align}
  F_\theta &= -\frac{1}{r}\ \frac{\partial u }{\partial \theta} = -J_2\ \frac{1}{r^4} 3\ \cos\theta\ \sin\theta \\
       F_r &= -\frac{\partial u }{\partial r} = J_2\ \frac{1}{r^4} \frac{3}{2}\ \left(3\sin^2\theta\ -\ 1\right)
\end{align}</math>|{{EquationRef|12}}}}

In the rectangular coordinate system (''x, y, z'') with unit vectors (''x̂ ŷ ẑ'') the force components are:
{{NumBlk|:|<math>\begin{align}
  F_x &= -\frac{\partial u }{\partial x} = J_2\ \frac{x}{r^7} \left(6z^2 - \frac{3}{2} (x^2 + y^2)\right) \\
  F_y &= -\frac{\partial u }{\partial y} = J_2\ \frac{y}{r^7} \left(6z^2 - \frac{3}{2} (x^2 + y^2)\right) \\
  F_z &= -\frac{\partial u }{\partial z} = J_2\ \frac{z}{r^7} \left(3z^2 - \frac{9}{2} (x^2 + y^2)\right)
\end{align}</math>|{{EquationRef|13}}}}

The components of the force corresponding to the "''J''<sub>3</sub> term"

:<math>u = \frac{J_3 P^0_3(\sin\theta) }{r^4} = J_3 \frac{1}{r^4} \frac{1}{2} \sin\theta \left(5\sin^2\theta - 3\right) = J_3 \frac{1}{r^7} \frac{1}{2} z \left(5 z^2 - 3 r^2\right)</math>

are
{{NumBlk|:|<math>\begin{align}
  F_\theta &= -\frac{1}{r} \frac{\partial u }{\partial \theta} = -J_3 \frac{1}{r^5} \frac{3}{2} \cos\theta \left(5 \sin^2\theta - 1\right) \\
       F_r &= -\frac{\partial u }{\partial r} = J_3 \frac{1}{r^5} 2 \sin\theta \left(5\sin^2\theta - 3\right)
\end{align}</math>|{{EquationRef|14}}}}

and

{{NumBlk|:|<math>\begin{align}
  F_x &= -\frac{\partial u }{\partial x} = J_3 \frac{x z}{r^9} \left(10 z^2 - \frac{15}{2} (x^2 + y^2)\right) \\
  F_y &= -\frac{\partial u }{\partial y} = J_3 \frac{y z}{r^9} \left(10 z^2 - \frac{15}{2} (x^2 + y^2)\right) \\
  F_z &= -\frac{\partial u }{\partial z} = J_3 \frac{1}{r^9} \left(4 z^2\ \left( z^2 - 3 (x^2 + y^2)\right) + \frac{3}{2} (x^2 + y^2)^2\right)
\end{align}</math>|{{EquationRef|15}}}}

The exact numerical values for the coefficients deviate (somewhat) between different Earth models but for the lowest coefficients they all agree almost exactly.

For the JGM-3 model (see below) the values are:

: μ = 398600.440 km<sup>3</sup>⋅s<sup>−2</sup>
: ''J''<sub>2</sub> = 1.75553 × 10<sup>10</sup> km<sup>5</sup>⋅s<sup>−2</sup>
: ''J''<sub>3</sub> = −2.61913 × 10<sup>11</sup> km<sup>6</sup>⋅s<sup>−2</sup>

For example, at a radius of 6600&nbsp;km (about 200&nbsp;km above Earth's surface) ''J''<sub>3</sub>/(''J''<sub>2</sub>''r'') is about 0.002; i.e., the correction to the "''J''<sub>2</sub> force" from the "''J''<sub>3</sub> term" is in the order of 2 permille. The negative value of ''J''<sub>3</sub> implies that for a point mass in Earth's equatorial plane the gravitational force is tilted slightly towards the south due to the lack of symmetry for the mass distribution of Earth's "north–south".

===Derivation===
{{Collapse top|title=A compact derivation of the spherical harmonics used to model Earth's gravitational field.}}

The spherical harmonics are derived from the approach of looking for harmonic functions of the form

{{NumBlk|:|<math>
  \phi = R(r)\ \Theta(\theta)\ \Phi(\varphi)
</math>|{{EquationRef|16}}}}

where (''r'', θ, φ) are the spherical coordinates defined by the equations ({{EquationNote|8}}). By straightforward calculations one gets that for any function ''f''

{{NumBlk|:|<math>
  \frac{\partial^2 f}{\partial x^2}\ +\ \frac{\partial^2 f}{\partial y^2}\ +\ \frac{\partial^2 f}{\partial z^2}\ =\ {1 \over r^2}{\partial \over \partial r}\left(r^2 {\partial f \over \partial r}\right) 
  + {1 \over r^2\cos\theta}{\partial \over \partial \theta}\left(\cos\theta {\partial f \over \partial \theta}\right) 
  + {1 \over r^2\cos^2\theta}{\partial^2 f \over \partial \varphi^2}
</math>|{{EquationRef|17}}}}

Introducing the expression ({{EquationNote|16}}) in ({{EquationNote|17}}) one gets that

{{NumBlk|:|<math>
\frac{r^2}{\phi}\left(\frac{\partial^2 \phi}{\partial x^2}+\frac{\partial^2 \phi}{\partial y^2} + \frac{\partial^2 \phi}{\partial z^2}\right)\ =
\frac{1}{R}\frac{d}{dr}\left(r^2\frac{dR}{dr}\right) + \frac{1}{\Theta\cos\theta}\frac{d}{d\theta}\left(\cos\theta \frac{d\Theta}{d\theta}\right) + \frac{1}{\Phi\cos^2\theta}\frac{d^2\Phi}{d\varphi^2} 
</math>|{{EquationRef|18}}}}

As the term 
:<math>\frac{1}{R}\frac{d}{dr}\left(r^2\frac{dR}{dr}\right)</math>
only depends on the variable <math>r</math> and the sum

:<math>\frac{1}{\Theta\cos\theta}\frac{d}{d\theta}\left(\cos\theta \frac{d\Theta}{d\theta}\right) + \frac{1}{\Phi\cos^2\theta}\frac{d^2\Phi}{d\varphi^2}</math>

only depends on the variables θ and φ. One gets that φ is harmonic if and only if 
{{NumBlk|:|<math>
\frac{1}{R}\frac{d}{dr}\left(r^2\frac{dR}{dr}\right) = \lambda
</math>|{{EquationRef|19}}}}

and
{{NumBlk|:|<math>
\frac{1}{\Theta\cos\theta}\frac{d}{d\theta}\left(\cos\theta \frac{d\Theta}{d\theta}\right) + \frac{1}{\Phi\cos^2\theta}\frac{d^2\Phi}{d\varphi^2} = -\lambda
</math>|{{EquationRef|20}}}}

for some constant <math>\lambda</math>.

From ({{EquationNote|20}}) then follows that

:<math>\frac{1}{\Theta}\ \cos\theta\ \frac{d}{d\theta}\left(\cos\theta \frac{d\Theta}{d\theta}\right)\ + \lambda\ \cos^2\theta\ +\ \frac{1}{\Phi}\frac{d^2\Phi}{d\varphi^2}\ =\ 0</math>

The first two terms only depend on the variable <math>\theta</math> and the third only on the variable <math>\varphi</math>.

From the definition of φ as a spherical coordinate it is clear that Φ(φ) must be periodic with the period 2π and one must therefore have that

{{NumBlk|:|<math>
\frac{1}{\Phi}\frac{d^2\Phi}{d\varphi^2} = -m^2
</math>|{{EquationRef|21}}}}

and
{{NumBlk|:|<math>
\frac{1}{\Theta}\ \cos\theta\ \frac{d}{d\theta}\left(\cos\theta \frac{d\Theta}{d\theta}\right)\ + \lambda\ \cos^2\theta = m^2
</math>|{{EquationRef|22}}}}

for some integer ''m'' as the family of solutions to ({{EquationNote|21}}) then are 
{{NumBlk|:|<math>
\Phi(\varphi)\ =a\ \cos m\varphi\ +\ b\ \sin m\varphi
</math>|{{EquationRef|23}}}}

With the variable substitution 
:<math>x=\sin \theta</math>

equation ({{EquationNote|22}}) takes the form

{{NumBlk|:|<math>
\frac{d}{dx}\left(\left(1 - x^2\right) \frac{d\Theta}{dx}\right) + \left(\lambda -\frac{m^2}{1 - x^2}\right)\Theta = 0
</math>|{{EquationRef|24}}}}

From ({{EquationNote|19}}) follows that in order to have a solution <math>\phi</math> with
:<math>R(r) = \frac{1}{r^{n+1}}</math>

one must have that

:<math>\lambda = n (n + 1)</math>

If ''P<sub>n</sub>''(''x'') is a solution to the differential equation
{{NumBlk|:|<math>
\frac{d}{dx}\left(\left(1 - x^2\right)\ \frac{dP_n}{dx}\right)\ +\ n(n + 1)\ P_n = 0
</math>|{{EquationRef|25}}}}

one therefore has that the potential corresponding to ''m'' = 0

:<math>\phi = \frac{1}{r^{n+1}}\ P_n(\sin\theta)</math>

'''which is rotationally symmetric around the z-axis''' is a harmonic function

If <math>P_{n}^{m}(x)</math> is a solution to the differential equation
{{NumBlk|:|<math>
\frac{d}{dx}\left(\left(1 - x^2\right)\ \frac{dP_{n}^{m}}{dx}\right)\ +\ \left(n(n + 1) - \frac{m^2}{1 - x^2} \right)\ P_{n}^{m}\ =\ 0
</math>|{{EquationRef|26}}}}

with ''m'' ≥ 1 one has the potential

{{NumBlk|:|<math>
\phi = \frac{1}{r^{n+1}}\ P_{n}^{m}(\sin\theta)\ (a\ \cos m\varphi\ +\ b\ \sin m\varphi)
</math>|{{EquationRef|27}}}}

where ''a'' and ''b'' are arbitrary constants is a harmonic function that depends on φ and therefore is '''not''' rotationally symmetric around the z-axis

The differential equation ({{EquationNote|25}}) is the Legendre differential equation for which the [[Legendre polynomials]] defined
{{NumBlk|:|<math>\begin{align}
  P_{0}(x) &= 1 \\
  P_{n}(x) &= \frac{1}{2^n n!}\ \frac{d^n\left(x^2 - 1\right)^n}{dx^n} \quad n \ge 1 \\
\end{align}</math>|{{EquationRef|28}}}}

are the solutions.

The arbitrary factor 1/(2<sup>n</sup>''n''!) is selected to make {{nowrap|''P<sub>n</sub>''(−1) {{=}} −1}} and {{nowrap|''P<sub>n</sub>''(1) {{=}} 1}} for odd ''n'' and {{nowrap|''P<sub>n</sub>''(−1) {{=}} ''P<sub>n</sub>''(1) {{=}} 1}} for even ''n''.

The first six Legendre polynomials are:

{{NumBlk|:|<math>\begin{align}
  P_{0}(x) &= 1 \\
  P_{1}(x) &= x \\
  P_{2}(x) &= \frac{1}{2} \left(3x^2 - 1\right) \\
  P_{3}(x) &= \frac{1}{2} \left(5x^3 - 3x\right) \\
  P_{4}(x) &= \frac{1}{8} \left(35x^4 - 30x^2 + 3\right) \\
  P_{5}(x) &= \frac{1}{8} \left(63x^5 - 70x^3 + 15x\right) \\
\end{align}</math>|{{EquationRef|29}}}}

The solutions to differential equation ({{EquationNote|26}}) are the associated [[Legendre polynomials|Legendre functions]]
{{NumBlk|:|<math>
  P_n^m(x) = \left(1 - x^2\right)^\frac{m}{2}\ \frac{d^m P_n}{dx^m} \quad 1 \le m \le n
</math>|{{EquationRef|30}}}}

One therefore has that

:<math>
  P_n^m (\sin\theta) = \cos^m \theta\ \frac{d^n P_n}{dx^n} (\sin\theta)
</math>
{{Collapse bottom}}

==Recursive algorithms used for the numerical propagation of spacecraft orbits==

Spacecraft orbits are computed by the [[Numerical methods for ordinary differential equations|numerical integration]] of the [[Physics:Equations of motion|equation of motion]]. For this the gravitational force, i.e. the [[Gradient|gradient]] of the potential, must be computed.  Efficient recursive algorithms have been designed to compute the gravitational force for any <math>N_z</math> and <math>N_t</math> (the max degree of zonal and tesseral terms) and such algorithms are used in standard orbit propagation software.

==Available models==

The earliest Earth models in general use by [[Organization:NASA|NASA]] and [[Organization:European Space Research Organisation|ESRO]]/[[Organization:European Space Agency|ESA]] were the "Goddard Earth Models" developed by [[Organization:Goddard Space Flight Center|Goddard Space Flight Center]] denoted "GEM-1", "GEM-2", "GEM-3", and so on. Later the "Joint Earth Gravity Models" denoted "JGM-1", "JGM-2", "JGM-3" developed by [[Organization:Goddard Space Flight Center|Goddard Space Flight Center]] in cooperation with universities and private companies became available. The newer models generally provided higher order terms than their precursors. The [[Earth:EGM96|EGM96]] uses ''N<sub>z</sub>'' = ''N<sub>t</sub>'' = 360 resulting in 130317 coefficients. An EGM2008 model is available as well.

For a normal Earth satellite requiring an orbit determination/prediction accuracy of a few meters the "JGM-3" truncated to ''N<sub>z</sub>'' = ''N<sub>t</sub>'' = 36 (1365 coefficients) is usually sufficient. Inaccuracies from the modeling of the air-drag and to a lesser extent the solar radiation pressure will exceed the inaccuracies caused by the gravitation modeling errors.

The dimensionless coefficients <math>\tilde{J_n} = -\frac{J_n}{\mu\ R^n}</math>, <math>\tilde{C_{n}^m} = -\frac{C_{n}^m}{\mu\ R^n}</math>, <math>\tilde{S_{n}^m} = -\frac{S_{n}^m}{\mu\ R^n}</math> for the first zonal and tesseral terms (using <math>R</math> = {{val|6378.1363|u=km}} and <math>\mu</math> = {{val|398600.4415|u=km<sup>3</sup>/s<sup>2</sup>}}) of the JGM-3 model are

{| class="wikitable"
|+ style="text-align: left;" | Zonal coefficients
|-
! scope="col" | n
!
|-
! scope="row" | 2  
| {{val|−0.1082635854E−2}}
|-
! scope="row" | 3
| {{val|0.2532435346E−5}}
|-
! scope="row" | 4
| {{val|0.1619331205E−5}}
|-
! scope="row" | 5
| {{val|0.2277161016E−6}}
|-
! scope="row" | 6
| {{val|−0.5396484906E−6}}
|-
! scope="row" | 7
| {{val|0.3513684422E−6}}
|-
! scope="row" | 8
| {{val|0.2025187152E−6}}
|}
{| class="wikitable"
|+ style="text-align: left;" |  Tesseral coefficients
|-
! scope="col" | ''n''
! scope="col" | ''m''
! scope="col" | ''C''
! scope="col" | ''S''
|-
! rowspan=2 scope="row" | 2
! scope="row" | 1 
| {{val|−0.3504890360E−9}} || {{val|0.1635406077E−8}}
|-
! scope="row" | 2 
| {{val|0.1574536043E−5}} || {{val|−0.9038680729E−6}}
|-
! rowspan=3 scope="row" | 3
! scope="row" | 1 
| {{val|0.2192798802E−5}} || {{val|0.2680118938E−6}}
|-
! scope="row" | 2 
| {{val|0.3090160446E−6}} || {{val|−0.2114023978E−6}}
|-
! scope="row" | 3 
| {{val|0.1005588574E−6}} || {{val|0.1972013239E−6}}
|-
! rowspan=4 scope="row" | 4
! scope="row" | 1 
| {{val|−0.5087253036E−6}} || {{val|−0.4494599352E−6}}
|-
! scope="row" | 2 
| {{val|0.7841223074E−7}} || {{val|0.1481554569E−6}}
|-
! scope="row" | 3 
| {{val|0.5921574319E−7}} || {{val|−0.1201129183E−7}} 
|-
! scope="row" | 4 
| {{val|−0.3982395740E−8}} || {{val|0.6525605810E−8}}
|}

According to JGM-3 one therefore has that ''J''{{sub|2}} = {{val|0.1082635854e−2}} × 6378.1363{{sup|2}} × {{val|398600.4415|u=km<sup>5</sup>/s<sup>2</sup>}} = {{val|1.75553e10|u=km<sup>5</sup>/s<sup>2</sup>}} and
''J''{{sup|3}} = {{val|-0.2532435346e−5}} × 6378.1363{{sup|3}} × {{val|398600.4415|u=km<sup>6</sup>/s<sup>2</sup>}} = {{val|−2.61913e11|u=km<sup>6</sup>/s<sup>2</sup>}}.

==Further reading==
* El'Yasberg ''Theory of flight of artificial earth satellites'', Israel program for Scientific Translations (1967)
* Lerch, F.J., Wagner, C.A., Smith, D.E., Sandson, M.L., Brownd, J.E., Richardson, J.A.,"Gravitational Field Models for the Earth (GEM1&2)", Report X55372146, Goddard Space Flight Center, Greenbelt/Maryland, 1972
* Lerch, F.J., Wagner, C.A., Putney, M.L., Sandson, M.L., Brownd, J.E., Richardson, J.A., Taylor, W.A., "Gravitational Field Models GEM3 and 4", Report X59272476, Goddard Space Flight Center, Greenbelt/Maryland, 1972
* Lerch, F.J., Wagner, C.A., Richardson, J.A., Brownd, J.E., "Goddard Earth Models (5 and 6)", Report X92174145, Goddard Space Flight Center, Greenbelt/Maryland, 1974
*  Lerch, F.J., Wagner, C.A., Klosko, S.M., Belott, R.P., Laubscher, R.E., Raylor, W.A., "Gravity Model Improvement Using Geos3 Altimetry (GEM10A and 10B)", 1978 Spring Annual Meeting of the American Geophysical Union, Miami, 1978
* Lerch, F.J., Klosko, S.M., Laubscher, R.E., Wagner, C.A., "Gravity Model Improvement Using Geos3 (GEM9 and 10)", Journal of Geophysical Research, Vol. 84, B8, p.&nbsp;3897-3916, 1979
* Lerch, F.J., Putney, B.H., Wagner, C.A., Klosko, S.M. ,"Goddard earth models for oceanographic applications (GEM 10B and 10C)", Marine-Geodesy, 5(2), p.&nbsp;145-187, 1981
* Lerch, F.J., Klosko, S.M., Patel, G.B., "A Refined Gravity Model from Lageos (GEML2)", 'NASA Technical Memorandum 84986, Goddard Space Flight Center, Greenbelt/Maryland, 1983
* Lerch, F.J., Nerem, R.S., Putney, B.H., Felsentreger, T.L., Sanchez, B.V., Klosko, S.M., Patel, G.B., Williamson, R.G., Chinn, D.S., Chan, J.C., Rachlin, K.E., Chandler, N.L., McCarthy, J.J., Marshall, J.A., Luthcke, S.B., Pavlis, D.W., Robbins, J.W., Kapoor, S., Pavlis, E.C., " Geopotential Models of the Earth from Satellite Tracking, Altimeter and Surface Gravity Observations: GEMT3 and GEMT3S", NASA Technical Memorandum 104555, Goddard Space Flight Center, Greenbelt/Maryland, 1992
* Lerch, F.J., Nerem, R.S., Putney, B.H., Felsentreger, T.L., Sanchez, B.V., Marshall, J.A., Klosko, S.M., Patel, G.B., Williamson, R.G., Chinn, D.S., Chan, J.C., Rachlin, K.E., Chandler, N.L., McCarthy, J.J., Luthcke, S.B., Pavlis, N.K., Pavlis, D.E., Robbins, J.W., Kapoor, S., Pavlis, E.C., "A Geopotential Model from Satellite Tracking, Altimeter and Surface Gravity Data: GEMT3", Journal of Geophysical Research, Vol. 99, No. B2, p.&nbsp;2815-2839, 1994
* Nerem, R.S., Lerch, F.J., Marshall, J.A., Pavlis, E.C., Putney, B.H., Tapley, B.D., Eanses, R.J., Ries, J.C., Schutz, B.E., Shum, C.K., Watkins, M.M., Klosko, S.M., Chan, J.C., Luthcke, S.B., Patel, G.B., Pavlis, N.K., Williamson, R.G., Rapp, R.H., Biancale, R., Nouel, F., "Gravity Model Developments for Topex/Poseidon: Joint Gravity Models 1 and 2", Journal of Geophysical Research, Vol. 99, No. C12, p.&nbsp;24421-24447, 1994a
* Tapley, B. D.  M. M. Watkins, J. C. Ries, G. W. Davis, R. J. Eanes, S. R. Poole, H. J. Rim, B. E. Schutz, C. K. Shum, R. S. Nerem, F. J. Lerch, J. A. Marshall, S. M. Klosko, N. K. Pavlis, and R. G. Williamson, “The Joint Gravity Model 3,” J. Geophys. Res., vol. 101, No. B12, ,December 1996

== External links ==

* http://cddis.nasa.gov/lw13/docs/papers/sci_lemoine_1m.pdf
* http://geodesy.geology.ohio-state.edu/course/refpapers/Tapley_JGR_JGM3_96.pdf

[[Category:Spaceflight concepts]]
[[Category:Gravity]]
[[Category:Earth orbits]]

{{Sourceattribution|Geopotential model}}