High-Precision Orbit Propagator

High-precision satellite propagation using Runga-Kutta methods for differential equation solving

Force models are mostly pulled from Montenbruck & Gill: https://link.springer.com/book/10.1007/978-3-642-58351-3

The propagator can also compute the state transition matrix, meaning position and velocity covariances can be propagated as well.

The default propagator uses a Runga-Kutta 9(8) integrator with coefficient computed by Verner: https://www.sfu.ca/~jverner/

This works much better than lower-order Runga-Kutta solvers such as Dormund-Prince, and I don’t know why it isn’t more popular in numerical packages

This module includes a function to propagate a position and time directly, and a convenience “satstate” object that represents satellite position, velocity, and optionally covariance and can propagate itself to different times

Forces included in the propagator:

  1. Earth gravity with higher-order zonal terms

  2. Sun, Moon gravity

  3. Radiation pressure

  4. Atmospheric drag: NRL-MISE 2000 density model, with option to include space weather effects (can be large)

class satkit.satprop.PropResult

PropResult is a dictionary containing the results of a high-precision orbit propagation that is returned by the “propagate” function

Phi

6x6 State transition matrix corresponding to each time. Output is Nx6x6 numpy matrix, where N is the lenght of the output “time” list. Not included if output_phi kwarg is set to false (the default)

pos

GCRF position in meters at output times. Output is Nx3 numpy matrix, where N is the number of times

Type:

npt.ArrayLike[float]

stats

Python dictionary with statistics for the propagation. This includes:

  • num_eval: Number of function evaluations of the force model required to get solution with desired accuracy

  • accepted_steps: Accepted steps in the adpative Runga-Kutta solver

  • rejected_steps: Rejected steps in the adaptive Runga-Kutta solver

Type:

(PropStats)

time

List of satkit.time objects at which time is computed

Type:

npt.ArrayLike[satkit.time]

vel

GCRF velocity in meters per second at output times. Output is Nx3 numpy matrix, where N is the number of times

Type:

npt.ArrayLike[float]

class satkit.satprop.PropStats

PropStats is a dictionary containing statistics for the propagation.

accepted_steps

Accepted steps in the adaptive Runga-Kutta solver

Type:

int

num_eval

Number of function evaluations of the force model computed during propagation function call

Type:

int

rejected_steps

Rejected steps in the adaptive Runga-Kutta solver

Type:

int

satkit.satprop.propagate(pos: numpy.typing.ArrayLike[float], vel: numpy.typing.ArrayLike[float], start: satkit.time, **kwargs)

High-precision orbit propagator

Notes

  • Propagator uses advanced Runga-Kutta integrators and includes the following forces:

    • Earth gravity with higher-order zonal terms

    • Sun, Moon gravity

    • Radiation pressure

    • Atmospheric drag: NRL-MISE 2000 density model, with option to include space weather effects (can be large)

  • Stop time must be set by keyword argument, either explicitely or by duration

  • Solid Earth tides are not (yet) included in the model

Propagate statellite ephemeris (position, velocity in gcrs & time) to new time and output new position and velocity via Runge-Kutta integration.

Inputs and outputs are all in the Geocentric Celestial Reference Frame (GCRF)

Parameters:

  • pos (npt.ArrayLike[float]) – 3-element numpy array representing satellite GCRF position in meters

  • vel (npt.ArrayLike[float]) – 3-element numpy array representing satellite GCRF velocity in m/s

  • tm (satkit.time) – satkit.time object representing instant at which satellite is at “pos” & “vel”

Keyword Arguments:

  • stop_time (satkit.time) – satkit.time object representing instant at which new position and velocity will be computed

  • duration_secs (float) – duration in seconds from “tm” for at which new position and velocity will be computed.

  • duration_days (float) – duration in days from “tm” at which new position and velocity will be computed.

  • duration (satkit.duration) – duration from “tm” at which new position & velocity will be computed.

  • dt_secs (float) – Interval in seconds between “starttime” and “stoptime” at which solution will also be computed

  • dt_days (float) – Interval in days between “starttime” and “stoptime” at which solution will also be computed

  • dt (satkit.duration) – Interval over which new position & velocity will be computed

  • output_phi (bool) – Output 6x6 state transition matrix between “starttime” and “stoptime” (and at intervals, if specified)

  • propsettings (propsettings) – “propsettings” object with input settings for the propagation. if left out, default will be used.

  • satproperties (satproperties_static) – “SatPropertiesStatic” object with drag and radiation pressure succeptibility of satellite. If left out, drag and radiation pressure are neglected

Returns:

Python dictionary with the following elements:

  • ”time”: list of satkit.time objects at which solution is computed

  • ”pos”: GCRF position in meters at “time”. Output is a Nx3 numpy matrix, where N is the length of the output “time” list

  • ”vel”: GCRF velocity in meters / second at “time”. Output is a Nx3 numpy matrix, where N is the length of the output “time” list

  • ”Phi”: 6x6 State transition matrix corresponding to each time. Output is Nx6x6 numpy matrix, where N is the lenght of the output “time” list. Not included if output_phi kwarg is set to false (the default)

  • ”stats”: Python dictionary with statistics for the propagation. This includes:

    • ”num_eval”: Number of function evaluations of the force model required to get solution with desired accuracy

    • ”accepted_steps”: Accepted steps in the adpative Runga-Kutta solver

    • ”rejected_steps”: Rejected steps in the adaptive Runga-Kutta solver

Return type:

(PropResult)

class satkit.satprop.propsettings

This class contains settings used in the high-precision orbit propgator part of the “satkit” python toolbox

property abs_error

Maxmum absolute value of error for any element in propagated state following ODE integration

Default: 1e-8

property gravity_order

Earth gravity order to use in ODE integration

Default: 4

property rel_error

Maximum relative error of any element in propagated state following ODE integration

Default: 1e-8

property use_jplephem

Use high-precision but computationally expensive JPL ephemerides for sun and mun when computing their gravitational force

property use_spaceweather

Use space weather data when computing atmospheric density for drag forces

Default: true

Notes

  • Space weather data can have a large effect on the density of the atmosphere

  • This can be important for accurate drag force calculations

  • Space weather data is updated every 3 hours. Most-recent data can be downloaded with satkit.utils.update_datafiles()

class satkit.satprop.satproperties_static(*args, **kwargs)

Satellite properties relevant for drag and radiation pressure

This class lets the satellite radiation pressure and drag paramters be set to static values for duration of propagation

class satkit.satprop.satstate(time: satkit.time, pos: numpy.typing.ArrayLike[numpy.float64], vel: numpy.typing.ArrayLike[numpy.float64], cov: numpy.typing.ArrayLike[numpy.float64] | None = None)

A convenience class representing a satellite position and velocity, and optionally 6x6 position/velocity covariance at a particular instant in time

This class can be used to propagate the position, velocity, and optional covariance to different points in time.

property cov

6x6 state covariance matrix in GCRF frame

Returns:

6x6 numpy array representing state covariance in GCRF frame or None if not set

Return type:

npt.ArrayLike[np.float64] | None

property pos_gcrf

state position in meters in GCRF frame

Returns:

3-element numpy array representing position in meters in GCRF frame

Return type:

npt.ArrayLike[np.float64]

propagate(time: satkit.time | satkit.duration, propsettings=None)

Propagate this state to a new time, specified by the “time” input, updating the position, the velocity, and the covariance if set

Parameters:

time (satkit.time|satkit.duration) – Time or duration from current time to which to propagate the state

Keyword Arguments:

propsettings – satkit.satprop.propsettings object describing settings to use in the propagation. If omitted, default is used

Returns:

New satellite state object representing the state at the new time

Return type:

satstate

property qgcrf2lvlh

Quaternion that rotates from the GCRF to the LVLH frame for the current state

Returns:

Quaternion that rotates from the GCRF to the LVLH frame for the current state

Return type:

satkit.quaternion

property time

Return time of this satellite state

Returns:

Time instant of this state

Return type:

satkit.time

property vel_gcrf

Return this state velocity in meters / second in GCRF

Returns:

3-element numpy array representing velocity in meters / second in GCRF frame

Return type:

npt.ArrayLike[np.float64]