Skip to main content
Log in

First-order analytic propagation of satellites in the exponential atmosphere of an oblate planet

  • Original Article
  • Published:
Celestial Mechanics and Dynamical Astronomy Aims and scope Submit manuscript

Abstract

The paper offers the fully analytic solution to the motion of a satellite orbiting under the influence of the two major perturbations, due to the oblateness and the atmospheric drag. The solution is presented in a time-explicit form, and takes into account an exponential distribution of the atmospheric density, an assumption that is reasonably close to reality. The approach involves two essential steps. The first one concerns a new approximate mathematical model that admits a closed-form solution with respect to a set of new variables. The second step is the determination of an infinitesimal contact transformation that allows to navigate between the new and the original variables. This contact transformation is obtained in exact form, and afterwards a Taylor series approximation is proposed in order to make all the computations explicit. The aforementioned transformation accommodates both perturbations, improving the accuracy of the orbit predictions by one order of magnitude with respect to the case when the atmospheric drag is absent from the transformation. Numerical simulations are performed for a low Earth orbit starting at an altitude of 350 km, and they show that the incorporation of drag terms into the contact transformation generates an error reduction by a factor of 7 in the position vector. The proposed method aims at improving the accuracy of analytic orbit propagation and transforming it into a viable alternative to the computationally intensive numerical methods.

This is a preview of subscription content, log in via an institution to check access.

Access this article

Price excludes VAT (USA)
Tax calculation will be finalised during checkout.

Instant access to the full article PDF.

Institutional subscriptions

Fig. 1
Fig. 2

Similar content being viewed by others

Notes

  1. http://www.lockheedmartin.com/us/products/space-fence.html.

  2. Due to the presence of significant uncertainties in the atmosphere (influence of the solar activity, day/night density variations, the atmospheric bulge, winds), realistic long-term propagations in the presence of drag cannot be addressed with deterministic tools (see Dell’Elce and Kerschen 2014).

  3. The critical inclination is a resonance of the dynamical system (see Lara 2015b) and therefore cannot be removed.

  4. Deprit was not the first one to use this approach. Hori (1966) proposed a slightly different version a few years before.

  5. They are uniquely defined up to a translation on the real numbers axis that is an integer multiple of \(2\pi \). But since only their trigonometric functions \(\sin \) and \(\cos \) are involved in computations, no restrictions need to be imposed.

  6. A transformation is said to be canonical if its Jacobian is a symplectic matrix, i.e. it obeys the condition in Eq. (19).

  7. Whittaker variables are singular for equatorial orbits, but this does not pose a problem for the case where only even zonal harmonics are taken into consideration, which is the case of the present paper.

  8. Elementary transformations between different sets of coordinates are also required at some steps of the algorithm, but these coordinates are present in the analytic expressions, making this requirement obvious,.

  9. These values were inspired by the prospective QARMAN satellite, a CubeSat mission developped by the von Karman Institute of Fluid Dynamics in Brussels and the University of Liège.

  10. The argument of perigee drift was omitted, since it is incorporated in the argument of latitude, depicted in Fig. 2c.

  11. Consequently, the orthogonal right-handed frame associated to the Poincaré canonical variables is found to be the same as for the equinoctial ones (see Battin 1999, Problem 10-10, p. 494).

References

  • Barrio, R., Palácian, J.: Lie Transforms for ordinary differential equations: taking advantage of the Hamiltonian form of terms of the perturbation. Int. J. Numer. Meth. Eng. 40(12), 2289–2300 (1997)

    Article  MathSciNet  MATH  Google Scholar 

  • Barrio, R., Palacián, J.: High-order averaging of eccentric artificial satellites perturbed by the earth’potential and air-drag terms. Proc. R. Soc. Lond. A Math. Phys. Eng. Sci. 459(2034), 1517–1534 (2003)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Battin, R.H.: An Introduction to the Mathematics and Methods of Astrodynamics. AIAA, Reston, VA (1999)

    Book  MATH  Google Scholar 

  • Brouwer, D.: Solution of the problem of artificial satellite theory without drag. Astron. J. 64, 378–396 (1959)

    Article  ADS  MathSciNet  Google Scholar 

  • Brouwer, D., Hori, G.I.: Theoretical evaluation of atmospheric drag effects in the motion of an artificial satellite. Astron. J. 66, 193–225 (1961)

    Article  ADS  MathSciNet  Google Scholar 

  • Celletti, A., Negrini, P.: Non-integrability of the problem of motion around an oblate planet. Celest. Mech. Dyn. Astron. 61, 253–260 (1995)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Cid, R., Lahulla, J.F.: Perturbaciones de corto periodo en el movimiento de un satélite artificial, en función de las variables de Hill. Publicaciones de la Revista de la Academia de Ciencias de Zaragoza 24, 159–165 (1969)

    Google Scholar 

  • Dell’Elce, L., Kerschen, G.: Probabilistic assessment of the lifetime of low-earth-orbit spacecraft: uncertainty characterization. J. Guid. Control Dyn. 1, 1–13 (2014)

    Google Scholar 

  • Deprit, A.: Canonical transformations depending on a small parameter. Celest. Mech. 1, 12–30 (1969)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Deprit, A.: The elimination of the parallax in satellite theory. Celest. Mech. 24, 111–153 (1981)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Garfinkel, B.: The orbit of a satellite of an oblate planet. Astron. J. 64, 353–366 (1959)

    Article  ADS  MathSciNet  Google Scholar 

  • Goldtsein, H., Poole, ChP, Safko., J.L.: Classical Mechanics, 3rd edn. Addison-Wesley, Boston (2002)

    Google Scholar 

  • Gurfil, P., Lara, M.: Satellite onboard orbit propagation using Deprit’s radial intermediary. Celest. Mech. Dyn. Astron. 120(2), 217–232 (2014)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Hilton, C.G., Kuhlman, J.R.: Mathematical Models for the Space Defense Center. Philco-Ford Corporation Publication U-3871, pp. 17–28 (1966)

  • Hori, G.: Theory of general perturbation with unspecified canonical variables. Publ. Astron. Soc. Jpn. 18, 287–296 (1966)

    ADS  Google Scholar 

  • Irigoyen, M., Simo, C.: Non integrability of the \(J_2\) problem. Celest. Mech. Dyn. Astron. 55, 281–287 (1993)

    Article  ADS  MATH  Google Scholar 

  • Lane, M.H., Cranford, K.H.: American Institute of Aeronautics, Astronautics, and American Astronautical Society. In: An Improved Analytical Drag Theory for the Artificial Satellite Problem. AIAA paper. American Institute of Aeronautics and Astronautics, (1969). https://books.google.be/books?id=dxeItgAACAAJ

  • Lara, M.: Efficient formulation of the periodic corrections in Brouwers gravity solution. Math. Probl. Eng. 2015, 9 (2015)

    Article  MathSciNet  Google Scholar 

  • Lara, M.: Leo intermediary propagation as a feasible alternative to brouwer’s gravity solution. Adv. Space Res. 56(3), 367–376 (2015b). (Advances in Asteroid and Space Debris Science and Technology—Part 1)

    Article  ADS  Google Scholar 

  • Lara, M., Gurfil, P.: Integrable approximation of \(J_2\)-perturbed relative orbits. Celest. Mech. Dyn. Astron. 114(3), 229–254 (2012)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Lyddane, R.H.: Small eccentricities or inclinations in the Brouwer theory of the artificial satellite. Astron. J. 68, 555–558 (1963)

    Article  ADS  MathSciNet  Google Scholar 

  • Martinusi, V., Dell’Elce, L., Kerschen, G.: Analytic propagation of near-circular satellite orbits in the atmosphere of an oblate planet. Celest. Mech. Dyn. Astron. 123(1), 85–103 (2015)

    Article  ADS  MathSciNet  MATH  Google Scholar 

  • Steinberg, S.: Lie series and nonlinear ordinary differential equations. J. Math. Anal. Appl. 101(1), 39–63 (1984)

    Article  MathSciNet  MATH  Google Scholar 

  • Sterne, T.E.: The gravitational orbit of a satellite of an oblate planet. Astron. J. 63, 28–40 (1958)

    Article  ADS  MathSciNet  Google Scholar 

  • Vallado, D.A.: Fundamentals of Astrodynamics and Applications. Microcosm, Cleveland (2001)

    MATH  Google Scholar 

  • Vallado, D.A., Crawford, P., Hujsak, R., Kelso, T.S.: Revisiting Spacetrack Report #3: Rev 2. In: AIAA Astrodynamics Specialists Conference and Exhibit, August 2006 (2006)

  • Vinti, J.P.: Theory of the orbit of an artificial satellite with use of spheroidal coordinates. Astron. J. 65, 353–354 (1960)

    Article  ADS  Google Scholar 

  • von Zeipel, H.: Recherches sur le mouvement des petites planètes. Almqvist & Wiksells, Stockholm (1916)

    MATH  Google Scholar 

Download references

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Vladimir Martinusi.

Appendix

Appendix

1.1 Lie series approximation for \(\overline{L}\)

Although the variable \(\overline{L}\) has an explicit expression, depicted in Eq. (10a), the numerical implementation of the complex error function, as well as its inverse, rise the problem of dealing with differences of very large numbers. It would be much more convenient to have an accurate approximation for \(\overline{L}\), and this is obtained by writing the solution to Eq. (8d) with the help of Lie series. For a complete description of the expansion mechanism, see Steinberg (1984).

Denote:

$$\begin{aligned} \lambda =\sqrt{\mu \alpha };\quad \beta =\frac{\lambda }{\overline{L}_{0}} \end{aligned}$$

The IVP to solve is rewritten as:

$$\begin{aligned} \dot{\overline{L}}=-\mu C_{0}\exp \left( \frac{\overline{L}_{0}-\mu r_{0}}{ \mu \alpha }\right) \exp \left( -\frac{\overline{L}^{2}-\overline{L}_{0}^{2} }{\lambda ^{2}}\right) ,\quad \overline{L}\left( t_{0}\right) =\overline{L}_{0} \end{aligned}$$

and the change of variable:

$$\begin{aligned} z=\frac{\overline{L}^{2}-\overline{L}_{0}^{2}}{\lambda ^{2}} \end{aligned}$$
(49)

is performed. Denote

$$\begin{aligned} \xi _{0}=\overline{L}_{0}\frac{2\mu }{\lambda ^{2}}C_{0}\exp \left( \frac{ \overline{L}_{0}-\mu r_{0}}{\mu \alpha }\right) \end{aligned}$$

The IVP in z is:

$$\begin{aligned} \dot{z}=-\xi _{0}\sqrt{z\beta ^{2}+1}\exp \left( -z\right) ,~z\left( t_{0}\right) =0 \end{aligned}$$
(50)

Consider D to be the differential operator, defined for any arbitrary function \(\phi \):

$$\begin{aligned} D\phi =-\xi _{0}\sqrt{z\beta ^{2}+1}\exp \left( -z\right) \frac{\partial \phi }{\partial z} \end{aligned}$$
(51)

and define its formal powers as:

$$\begin{aligned} D^{k+1}\phi =D\left( D^{k}\phi \right) ,\;k\ge 1 \end{aligned}$$

Consider the formal exponential operator:

$$\begin{aligned} e^{\left( t-t_{0}\right) D}\phi =\phi +\sum _{k=1}^{\infty }\frac{\left( t-t_{0}\right) ^{k}}{k!}\left( D^{k}\phi \right) \end{aligned}$$
(52)

Then the solution to IVP (50) is expressed as:

$$\begin{aligned} z\left( t\right) =\left. e^{\left( t-t_{0}\right) D}z\right| _{z=0}, \end{aligned}$$
(53)

and explicitely:

$$\begin{aligned} z\left( t\right) =\sum _{k=1}^{\infty }\frac{\left( t-t_{0}\right) ^{k}}{k!} \left( \left. D^{k}z\right| _{z=0}\right) \end{aligned}$$
(54)

By taking into account the expression in Eq. (51) of the operator D,  the formal powers \(\left. D^{k}z\right| _{z=0}\) are expressed as:

$$\begin{aligned} \left. D^{1}z\right| _{z=0}= & {} -\xi _{0} \\ \left. D^{2}z\right| _{z=0}= & {} \xi _{0}^{2}\left( \frac{1}{2}\beta ^{2}-1\right) \\ \left. D^{3}z\right| _{z=0}= & {} \xi _{0}^{3}\left( 2\beta ^{2}-2\right) \\ \left. D^{4}z\right| _{z=0}= & {} \xi _{0}^{4}\left( -\beta ^{4}+9\beta ^{2}-6\right) \\ \left. D^{5}z\right| _{z=0}= & {} \xi _{0}^{5}\left( -13\beta ^{4}+48\beta ^{2}-24\right) \\ \left. D^{6}z\right| _{z=0}= & {} \xi _{0}^{6}\left( \frac{13}{2}\beta ^{6}-137\beta ^{4}+300\beta ^{2}-120\right) \\ \left. D^{7}z\right| _{z=0}= & {} \xi _{0}^{7}\left( 176\beta ^{6}-1422\beta ^{4}+2160\beta ^{2}-720\right) \\ \left. D^{8}z\right| _{z=0}= & {} \xi _{0}^{8}\left( -88\beta ^{8}+3365\beta ^{6}-15354\beta ^{4}+17640\beta ^{2}-5040\right) \\ \left. D^{9}z\right| _{z=0}= & {} \xi _{0}^{9}\left( -4069\beta ^{8}+57628\beta ^{6}-175752\beta ^{4}+161280\beta ^{2}-40320\right) \end{aligned}$$

Then \(\overline{L}=\overline{L}\left( t\right) \) is evaluated based on Eq. ( 49):

$$\begin{aligned} \overline{L}\left( t\right) =\overline{L}_{0}\sqrt{\beta ^{2}z\left( t\right) +1} \end{aligned}$$

1.2 Evaluation of \(I_{1,2,3}\) integrals

In the form they are written in Eq. (13), the integrals \(I_{1,2,3}\) are not suitable for an efficient series expansion approximation, given the considerable length of the interval on which the integration is performed (i.e., for a good approximation, a large number of terms are required in the truncated series expansion). For a rapid convergence, and also to avoid differences between very large numbers (since \(L=\sqrt{\mu \alpha }\) has the same order of magnitude of the angular momentum). It would therefore be ideal to transform them to integrals from 0 to some value z that has the smallest possible absolute value, so that the accuracy of the series expansion is maximum.

First, note that \(I_{1,2,3}\) may be rewritten as:

$$\begin{aligned} I_{k}=\exp \left( \frac{\overline{a}_{0}-r_{0}}{\alpha }\right) \int \limits _{ \overline{L}_{0}}^{\overline{L}}q_{k}\left( u\right) \exp \left( \frac{u^{2}- \overline{L}_{0}^{2}}{\mu \alpha }\right) \text {d}u,\quad k=1,2,3 \end{aligned}$$

where

$$\begin{aligned} q_{1}\left( u\right) =u^{-3},\quad q_{2}\left( u\right) =\left( u+\overline{G} _{0}-\overline{L}_{0}\right) ^{-3}u^{-4},\quad q_{3}\left( u\right) =\left( u+ \overline{G}_{0}-\overline{L}_{0}\right) ^{-4}u^{-3} \end{aligned}$$

To this end, denote

$$\begin{aligned} \lambda =\sqrt{\mu \alpha },\quad z=\frac{\overline{L}^{2}-\overline{L}_{0}^{2} }{\lambda ^{2}},\quad w_{k}\left( y\right) =\frac{q_{k}\left( \sqrt{\lambda ^{2}y+\overline{L}_{0}^{2}}\right) }{\sqrt{\lambda ^{2}y+\overline{L}_{0}^{2} }} \end{aligned}$$

and perform the change of variable:

$$\begin{aligned} y=\frac{u^{2}-\overline{L}_{0}^{2}}{\lambda ^{2}}\Rightarrow \text {d}u=\frac{ \lambda ^{2}}{2}\frac{\text {d}y}{\sqrt{\lambda ^{2}y+\overline{L}_{0}^{2}}} \end{aligned}$$

The integrals become:

$$\begin{aligned} I_{k}=\frac{\lambda ^{2}}{2}\exp \left( \frac{\overline{a}_{0}-r_{0}}{\alpha }\right) \int \limits _{0}^{z}w_{k}\left( y\right) \exp \left( y\right) \text {d}y \end{aligned}$$

By performing a Taylor series expansion for \(\exp \left( y\right) ,\) the integrals become:

$$\begin{aligned} I_{k}=\frac{\lambda ^{2}}{2}\exp \left( \frac{\overline{a}_{0}-r_{0}}{\alpha }\right) \sum \limits _{N=0}^{\infty }\left[ \frac{1}{N!}\int \limits _{0}^{z}y^{N}w_{k}\left( y\right) \text {d}y\right] \end{aligned}$$

and written separately:

$$\begin{aligned} I_{1}= & {} \frac{\lambda ^{2}}{2}\exp \left( \frac{\overline{a}_{0}-r_{0}}{ \alpha }\right) \sum \limits _{N=0}^{\infty }\left[ \frac{1}{N!} \int \limits _{0}^{z}\frac{y^{N}}{\left( \lambda ^{2}y+\overline{L} _{0}^{2}\right) ^{2}}\text {d}y\right] \\ I_{2}= & {} \frac{\lambda ^{2}}{2}\exp \left( \frac{\overline{a}_{0}-r_{0}}{ \alpha }\right) \sum \limits _{N=0}^{\infty }\left[ \frac{1}{N!} \int \limits _{0}^{z}\frac{y^{N}}{\left[ \sqrt{\lambda ^{2}y+\overline{L} _{0}^{2}}+\overline{G}_{0}-\overline{L}_{0}\right] ^{3}\left( \lambda ^{2}y+ \overline{L}_{0}^{2}\right) ^{\frac{5}{2}}}\text {d}y\right] \\ I_{3}= & {} \frac{\lambda ^{2}}{2}\exp \left( \frac{\overline{a}_{0}-r_{0}}{ \alpha }\right) \sum \limits _{N=0}^{\infty }\left[ \frac{1}{N!} \int \limits _{0}^{z}\frac{y^{N}}{\left[ \sqrt{\lambda ^{2}y+\overline{L} _{0}^{2}}+\overline{G}_{0}-\overline{L}_{0}\right] ^{4}\left( \lambda ^{2}y+ \overline{L}_{0}^{2}\right) ^{3}}\text {d}y\right] \end{aligned}$$

The evaluation of \(I_{1,2,3}\) has been simplified to definite integrals that may be expressed explicitly through elementary functions.

1.3 Jacobian for coordinate transformations

The Jacobian of the transformation from Delaunay to Whittaker elements is:

$$\begin{aligned} {\mathcal {J}}_{\mathcal {D}}^{{\mathcal {W}}}=\frac{\partial \left( r,\theta ,\nu ,R,\varTheta ,N\right) }{\partial \left( l,g,h,L,G,H\right) }=\left[ \begin{array}{cccccc} \dfrac{p\sigma }{\eta ^{3}} &{} 0 &{} 0 &{} \dfrac{r\left( \kappa ^{2}+\kappa -2e^{2}\right) }{Le^{2}} &{} \dfrac{\eta p\kappa }{Le^{2}} &{} 0 \\ &{} &{} &{} &{} &{} \\ \dfrac{\left( 1+\kappa \right) ^{2}}{\eta ^{3}} &{} 1 &{} 0 &{} \dfrac{\sigma \left( 2+\kappa \right) }{Le^{2}} &{} -\dfrac{\eta \sigma \left( 2+\kappa \right) }{Le^{2}} &{} 0 \\ &{} &{} &{} &{} &{} \\ 0 &{} 0 &{} 1 &{} 0 &{} 0 &{} 0 \\ &{} &{} &{} &{} &{} \\ \dfrac{L\kappa \left( 1+\kappa \right) ^{2}}{\eta ^{2}p} &{} 0 &{} 0 &{} \dfrac{ \eta \sigma }{pe^{2}}\left[ 1+2\kappa -\sigma ^{2}\right] &{} -\dfrac{\sigma \left( 1+\kappa \right) ^{2}}{pe^{2}} &{} 0 \\ &{} &{} &{} &{} &{} \\ 0 &{} 0 &{} 0 &{} 0 &{} 1 &{} 0 \\ &{} &{} &{} &{} &{} \\ 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 1 \end{array} \right] \end{aligned}$$

The Jacobian of the transformation from Poincaré to Delaunay variables is:

$$\begin{aligned} {\mathcal {J}}_{\mathcal {P}}^{\mathcal {D}}= & {} \frac{\partial \left( l,g,h,L,G,H\right) }{\partial \left( Q_{1},Q_{2},Q_{3},P_{1},P_{2},P_{3}\right) }\\= & {} \left[ \begin{array}{cccccc} 1 &{} \dfrac{P_{2}}{Q_{2}^{2}+P_{2}^{2}} &{} 0 &{} 0 &{} -\dfrac{Q_{2}}{ Q_{2}^{2}+P_{2}^{2}} &{} 0 \\ &{} &{} &{} &{} &{} \\ 0 &{} -\dfrac{P_{2}}{Q_{2}^{2}+P_{2}^{2}} &{} \dfrac{P_{3}}{Q_{3}^{2}+P_{3}^{2}} &{} 0 &{} \dfrac{Q_{2}}{Q_{2}^{2}+P_{2}^{2}} &{} -\dfrac{Q_{3}}{Q_{3}^{2}+P_{3}^{2} } \\ &{} &{} &{} &{} &{} \\ 0 &{} 0 &{} -\dfrac{P_{3}}{Q_{3}^{2}+P_{3}^{2}} &{} 0 &{} 0 &{} \dfrac{Q_{3}}{ Q_{3}^{2}+P_{3}^{2}} \\ &{} &{} &{} &{} &{} \\ 0 &{} 0 &{} 0 &{} 1 &{} 0 &{} 0 \\ &{} &{} &{} &{} &{} \\ 0 &{} -Q_{2} &{} 0 &{} 1 &{} -P_{2} &{} 0 \\ &{} &{} &{} &{} &{} \\ 0 &{} -Q_{2} &{} -Q_{3} &{} 1 &{} -P_{2} &{} -P_{3} \end{array} \right] \end{aligned}$$

The Jacobian of the transformation from Whittaker to Poincaré elements may now be explicitly derived as:

$$\begin{aligned} {\mathcal {J}}_{\mathcal {P}}^{{\mathcal {W}}}=\frac{\partial \left( r,\theta ,\nu ,R,\varTheta ,N\right) }{\partial \left( Q_{1},Q_{2},Q_{3},P_{1},P_{2},P_{3}\right) }={\mathcal {J}}_{\mathcal {D}}^{ {\mathcal {W}}}{\mathcal {J}}_{\mathcal {P}}^{\mathcal {D}} \end{aligned}$$

Since it involves only elementary computations, that are nevertheless leading to relatively long expressions, the full matrix form of \({\mathcal {J}} _{\mathcal {P}}^{{\mathcal {W}}}\) will not be displayed. In any case, it is obtained by the multiplication of two explicit matrices.

1.4 Transformation from Poincaré to Cartesian coordinates

The Poincaré canonical variables are mentioned and used in very few Astrodynamics/Orbital Mechanics textbooks, but they lack a comprehensive approach, and to our knowledge, a straightforward transformation from these variables to Cartesian is absent from the literature. The closest approach is found in Battin (1999), but it is made for the equinoctial variables, that are nonsingular but not canonical.

The transformation is determined here in the most direct way possible. It is known that the inertial position and velocity vectors are expressed in the perifocal frame \(\left\{ P\right\} \) as:

$$\begin{aligned} {\mathbf {r}}^{\left\{ P\right\} }=r\left[ \begin{array}{c} \cos f \\ \sin f \\ 0 \end{array} \right] ;\quad {\mathbf {v}}^{\left\{ P\right\} }=\frac{\mu }{G}\left[ \begin{array}{c} -\sin f \\ e+\cos f \\ 0 \end{array} \right] \end{aligned}$$

and their corresponding column matrices in the inertial frame originated in the attraction center are:

$$\begin{aligned} {\mathbf {r}}^{\left\{ IN\right\} }= & {} {\mathbf {R}}\left( {\mathbf {i}}_{3},h\right) {\mathbf {R}}\left( {\mathbf {i}}_{1},i\right) {\mathbf {R}}\left( {\mathbf {i}} _{3},g\right) {\mathbf {r}}^{\left\{ P\right\} } \\ {\mathbf {v}}^{\left\{ IN\right\} }= & {} {\mathbf {R}}\left( {\mathbf {i}}_{3},h\right) {\mathbf {R}}\left( {\mathbf {i}}_{1},i\right) {\mathbf {R}}\left( {\mathbf {i}} _{3},g\right) {\mathbf {v}}^{\left\{ P\right\} } \end{aligned}$$

where the orthogonal matrices \({\mathbf {R}}_{h},\) \({\mathbf {R}}_{i},\) \({\mathbf {R}} _{g}\) are:

$$\begin{aligned} {\mathbf {R}}\left( {\mathbf {i}}_{3},h\right)= & {} \left[ \begin{array}{ccc} \cos h &{} -\sin h &{} 0 \\ \sin h &{} \cos h &{} 0 \\ 0 &{} 0 &{} 1 \end{array} \right] ,\quad {\mathbf {R}}\left( {\mathbf {i}}_{1},i\right) =\left[ \begin{array}{ccc} 1 &{} 0 &{} 0 \\ 0 &{} \cos i &{} -\sin i \\ 0 &{} \sin i &{} \cos i \end{array} \right] ,\\ {\mathbf {R}}\left( {\mathbf {i}}_{3},g\right)= & {} \left[ \begin{array}{ccc} \cos g &{} -\sin g &{} 0 \\ \sin g &{} \cos g &{} 0 \\ 0 &{} 0 &{} 1 \end{array} \right] \end{aligned}$$

If a new fixed frame \(\left\{ \mathcal {P}\right\} \) is defined such that the position vector has the expression:

$$\begin{aligned} {\mathbf {r}}^{\left\{ \mathcal {P}\right\} }=r\left[ \begin{array}{c} \cos \left( f+g+h\right) \\ \sin \left( f+g+h\right) \\ 0 \end{array} \right] =r\left[ \begin{array}{c} \cos \left( Q_{1}+f-l\right) \\ \sin \left( Q_{1}+f-l\right) \\ 0 \end{array} \right] , \end{aligned}$$

it will follow (by elementary computations) that the velocity in the same frame \(\left\{ \mathcal {P}\right\} \) is:

$$\begin{aligned} {\mathbf {v}}^{\left\{ \mathcal {P}\right\} }=\frac{2\mu }{\sqrt{4P_{1}-2\left( Q_{2}^{2}+P_{2}^{2}\right) }}\left[ \begin{array}{c} \dfrac{Q_{2}}{2P_{1}}\sqrt{4P_{1}-\left( Q_{2}^{2}+P_{2}^{2}\right) }-\sin \left( Q_{1}+f-l\right) \\ \dfrac{P_{2}}{2P_{1}}\sqrt{4P_{1}-\left( Q_{2}^{2}+P_{2}^{2}\right) }+\cos \left( Q_{1}+f-l\right) \\ 0 \end{array} \right] \end{aligned}$$

Consequently, the inertial counterparts of \({\mathbf {r}}^{\left\{ \mathcal {P} \right\} },\) \({\mathbf {v}}^{\left\{ \mathcal {P}\right\} }\) are:

$$\begin{aligned} {\mathbf {r}}^{\left\{ IN\right\} }= & {} \mathbf {Ar}^{\left\{ \mathcal {P}\right\} } \\ {\mathbf {v}}^{\left\{ IN\right\} }= & {} \mathbf {Av}^{\left\{ \mathcal {P}\right\} } \end{aligned}$$

The matrix \(\mathbf {A}\) is found to beFootnote 11:

$$\begin{aligned} \mathbf {A}={\mathbf {R}}\left( {\mathbf {i}}_{3},h\right) {\mathbf {R}}\left( \mathbf { i}_{1},i\right) {\mathbf {R}}\left( {\mathbf {i}}_{3},-h\right) \end{aligned}$$

Define U such that:

$$\begin{aligned} U=\sqrt{4P_{1}-2\left( Q_{2}^{2}+P_{2}^{2}\right) -\left( Q_{3}^{2}+P_{3}^{2}\right) }=\sqrt{2G\left( 1+\cos i\right) } \end{aligned}$$

Then the explicit expression of matrix \(\mathbf {A,}\) in Poincaré variables, is:

$$\begin{aligned} \mathbf {A}=\frac{1}{U^{2}+Q_{3}^{2}+P_{3}^{2}}\left[ \begin{array}{ccc} U^{2}-Q_{3}^{2}+P_{3}^{2} &{} -2P_{3}Q_{3} &{} -2{{ UQ}}_{3} \\ -2P_{3}Q_{3} &{} U^{2}+Q_{3}^{2}-P_{3}^{2} &{} -2{{ UP}}_{3} \\ 2{{ UQ}}_{3} &{} 2{{ UP}}_{3} &{} U^{2}-Q_{3}^{2}-P_{3}^{2} \end{array} \right] \end{aligned}$$

Rights and permissions

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Martinusi, V., Dell’Elce, L. & Kerschen, G. First-order analytic propagation of satellites in the exponential atmosphere of an oblate planet. Celest Mech Dyn Astr 127, 451–476 (2017). https://doi.org/10.1007/s10569-016-9734-8

Download citation

  • Received:

  • Revised:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1007/s10569-016-9734-8

Keywords

Navigation