poliastro.core.propagation

Low level propagation algorithms

Package Contents

Functions

farnocchia(k, r0, v0, tof)

Propagates orbit using mean motion.

farnocchia_coe(k, p, ecc, inc, raan, argp, nu, tof)

func_twobody(t0, u_, k, ad, ad_kwargs)

Differential equation for the initial value two body problem.

vallado(k, r0, v0, tof, numiter)

Solves Kepler’s Equation by applying a Newton-Raphson method.

mikkola_coe(k, p, ecc, inc, raan, argp, nu, tof)

mikkola(k, r0, v0, tof, rtol=None)

Raw algorithm for Mikkola’s Kepler solver.

markley_coe(k, p, ecc, inc, raan, argp, nu, tof)

markley(k, r0, v0, tof)

Solves the kepler problem by a non iterative method. Relative error is

pimienta_coe(k, p, ecc, inc, raan, argp, nu, tof)

pimienta(k, r0, v0, tof)

Raw algorithm for Adonis’ Pimienta and John L. Crassidis 15th order

gooding_coe(k, p, ecc, inc, raan, argp, nu, tof, numiter=150, rtol=1e-08)

gooding(k, r0, v0, tof, numiter=150, rtol=1e-08)

Solves the Elliptic Kepler Equation with a cubic convergence and

danby_coe(k, p, ecc, inc, raan, argp, nu, tof, numiter=20, rtol=1e-08)

danby(k, r0, v0, tof, numiter=20, rtol=1e-08)

Kepler solver for both elliptic and parabolic orbits based on Danby’s

poliastro.core.propagation.farnocchia(k, r0, v0, tof)

Propagates orbit using mean motion.

This algorithm depends on the geometric shape of the orbit. For the case of the strong elliptic or strong hyperbolic orbits:

\[M = M_{0} + \frac{\mu^{2}}{h^{3}}\left ( 1 -e^{2}\right )^{\frac{3}{2}}t\]

New in version 0.9.0.

Parameters
  • k (float) – Standar Gravitational parameter

  • r0 (Quantity) – Initial position vector wrt attractor center.

  • v0 (Quantity) – Initial velocity vector.

  • tof (float) – Time of flight (s).

Note

This method takes initial \(\vec{r}, \vec{v}\), calculates classical orbit parameters, increases mean anomaly and performs inverse transformation to get final \(\vec{r}, \vec{v}\) The logic is based on formulae (4), (6) and (7) from http://dx.doi.org/10.1007/s10569-013-9476-9

poliastro.core.propagation.farnocchia_coe(k, p, ecc, inc, raan, argp, nu, tof)
poliastro.core.propagation.func_twobody(t0, u_, k, ad, ad_kwargs)

Differential equation for the initial value two body problem.

This function follows Cowell’s formulation.

Parameters
  • t0 (float) – Time.

  • u (ndarray) – Six component state vector [x, y, z, vx, vy, vz] (km, km/s).

  • k (float) – Standard gravitational parameter.

  • ad (function(t0, u, k)) – Non Keplerian acceleration (km/s2).

  • ad_kwargs (optional) – perturbation parameters passed to ad.

poliastro.core.propagation.vallado(k, r0, v0, tof, numiter)

Solves Kepler’s Equation by applying a Newton-Raphson method.

If the position of a body along its orbit wants to be computed for an specific time, it can be solved by terms of the Kepler’s Equation:

\[E = M + e\sin{E}\]

In this case, the equation is written in terms of the Universal Anomaly:

\[\sqrt{\mu}\Delta t = \frac{r_{o}v_{o}}{\sqrt{\mu}}\chi^{2}C(\alpha \chi^{2}) + (1 - \alpha r_{o})\chi^{3}S(\alpha \chi^{2}) + r_{0}\chi\]

This equation is solved for the universal anomaly by applying a Newton-Raphson numerical method. Once it is solved, the Lagrange coefficients are returned:

\[\begin{split}\begin{align} f &= 1 \frac{\chi^{2}}{r_{o}}C(\alpha \chi^{2}) \\ g &= \Delta t - \frac{1}{\sqrt{\mu}}\chi^{3}S(\alpha \chi^{2}) \\ \dot{f} &= \frac{\sqrt{\mu}}{rr_{o}}(\alpha \chi^{3}S(\alpha \chi^{2}) - \chi) \\ \dot{g} &= 1 - \frac{\chi^{2}}{r}C(\alpha \chi^{2}) \\ \end{align}\end{split}\]

Lagrange coefficients can be related then with the position and velocity vectors:

\[\begin{split}\begin{align} \vec{r} &= f\vec{r_{o}} + g\vec{v_{o}} \\ \vec{v} &= \dot{f}\vec{r_{o}} + \dot{g}\vec{v_{o}} \\ \end{align}\end{split}\]
Parameters
  • k (float) – Standard gravitational parameter

  • r0 (array) – Initial position vector

  • v0 (array) – Initial velocity vector

  • numiter (int) – Number of iterations

Returns

  • f (float) – First Lagrange coefficient

  • g (float) – Second Lagrange coefficient

  • fdot (float) – Derivative of the first coefficient

  • gdot (float) – Derivative of the second coefficient

Note

The theoretical procedure is explained in section 3.7 of Curtis in really deep detail. For analytical example, check in the same book for example 3.6.

poliastro.core.propagation.mikkola_coe(k, p, ecc, inc, raan, argp, nu, tof)
poliastro.core.propagation.mikkola(k, r0, v0, tof, rtol=None)

Raw algorithm for Mikkola’s Kepler solver.

Parameters
  • k (Quantity) – Standard gravitational parameter of the attractor.

  • r (Quantity) – Position vector.

  • v (Quantity) – Velocity vector.

  • tofs (Quantity) – Array of times to propagate.

  • rtol (float) – This method does not require of tolerance since it is non iterative.

Returns

  • rr (~astropy.units.Quantity) – Propagated position vectors.

  • vv (~astropy.units.Quantity)

Note

Original paper: https://doi.org/10.1007/BF01235850

poliastro.core.propagation.markley_coe(k, p, ecc, inc, raan, argp, nu, tof)
poliastro.core.propagation.markley(k, r0, v0, tof)

Solves the kepler problem by a non iterative method. Relative error is around 1e-18, only limited by machine double-precission errors.

Parameters
  • k (float) – Standar Gravitational parameter

  • r0 (array) – Initial position vector wrt attractor center.

  • v0 (array) – Initial velocity vector.

  • tof (float) – Time of flight.

Returns

  • rf (array) – Final position vector

  • vf (array) – Final velocity vector

Note

The following algorithm was taken from http://dx.doi.org/10.1007/BF00691917.

poliastro.core.propagation.pimienta_coe(k, p, ecc, inc, raan, argp, nu, tof)
poliastro.core.propagation.pimienta(k, r0, v0, tof)

Raw algorithm for Adonis’ Pimienta and John L. Crassidis 15th order polynomial Kepler solver.

Parameters
  • k (float) – Standar Gravitational parameter

  • r0 (array) – Initial position vector wrt attractor center.

  • v0 (array) – Initial velocity vector.

  • tof (float) – Time of flight.

Returns

  • rf (array) – Final position vector

  • vf (array) – Final velocity vector

Note

This algorithm was drived from the original paper: http://hdl.handle.net/10477/50522

poliastro.core.propagation.gooding_coe(k, p, ecc, inc, raan, argp, nu, tof, numiter=150, rtol=1e-08)
poliastro.core.propagation.gooding(k, r0, v0, tof, numiter=150, rtol=1e-08)

Solves the Elliptic Kepler Equation with a cubic convergence and accuracy better than 10e-12 rad is normally achieved. It is not valid for eccentricities equal or higher than 1.0.

Parameters
  • k (float) – Standard gravitational parameter of the attractor.

  • r (1x3 vector) – Position vector.

  • v (1x3 vector) – Velocity vector.

  • tof (float) – Time of flight.

  • rtol (float) – Relative error for accuracy of the method.

Returns

rr – Propagated position vectors. vv : 1x3 vector

Return type

1x3 vector

Note

Original paper for the algorithm: https://doi.org/10.1007/BF01238923

poliastro.core.propagation.danby_coe(k, p, ecc, inc, raan, argp, nu, tof, numiter=20, rtol=1e-08)
poliastro.core.propagation.danby(k, r0, v0, tof, numiter=20, rtol=1e-08)

Kepler solver for both elliptic and parabolic orbits based on Danby’s algorithm.

Parameters
  • k (float) – Standard gravitational parameter of the attractor.

  • r (1x3 vector) – Position vector.

  • v (1x3 vector) – Velocity vector.

  • tof (float) – Time of flight.

  • rtol (float) – Relative error for accuracy of the method.

Returns

  • rr (1x3 vector) – Propagated position vectors.

  • vv (1x3 vector)

Note

This algorithm was developed by Danby in his paper The solution of Kepler Equation with DOI: https://doi.org/10.1007/BF01686811