State Transition Matrix computation¶
Nyx uses dual number theory (i.e. automatic differentiation) to compute all state transition matrices (or Jacobian matrices) to machine precision without relying on finite differencing or errorprone manual derivations. All of the computations are handled by the hyperdual Rust package, which I coauthored.
The following is an extract of the AAS19716 conference paper called Application of Dual Number Theory to Statistical Orbital Determination, by C. Rabotin ^{1}.
First, the mathematical notions of dual numbers and hyperdual spaces is recalled, along with their errorfree autodifferentiation properties. Subsequently, hyperdual spaces are elegantly applied to statistical orbital determination problems, both for the computation of linearized dynamics, and for the sensitivity matrix of a Kalman filter. Further, a comparison of the execution speeds of the same statistical orbital determination problem using the hyperdual space formulation and the traditional manual derivation formulation is provided.
Please refer to Dual Numbers for a primer on dual number theory.
Applying hyperdual space to statistical orbital determination using Kalman filtering¶
The state transition matrix¶
The state transition matrix (STM, \(\boldsymbol \Phi\)) is a linearization procedure of a dynamical system. In the field of astrodynamics, it is an approximation of the dynamics over a short period of time. In the case of statistical orbital determination, the STM is computed around a reference point of an orbit and propagated forward in time until the next propagation time step or until the next spacecraft observation. The STM is computed via the partials matrix \(\boldsymbol A\) of the state \(\boldsymbol X\) at a time \(t\), as shown earlier.
If one only estimates the spacecraft state, then the partials matrix corresponds to the gradient of the orbital dynamics applied to the spacecraft. The state \(\boldsymbol X\) of a spacecraft may be defined as the vector of its position and velocity in a Cartesian frame. The equation below shows the relation between the partials matrix and the STM, where \(t\) and \(t_0\) correspond to two different times such that \(t>t_0\).
The \(\boldsymbol A\) matrix contains the partial derivatives of the accelerations (\(a_x,~a_y,~a_z\)) with respect to each of the components of the position of the spacecraft, cf. the next equation.
In the case of twobody dynamics, the \(\boldsymbol A\) matrix is relatively simple to compute. Deriving this partials matrix for higher fidelity dynamics or for a larger state is more complicated. In practice, this may lead orbital estimation software to ignore part of these dynamics. Navigators may then account for smaller perturbations by inflating the covariance matrix with stochastic noise.
Dual numbers, however, provide the partials matrix as part of the computation of the equations of motion (EOM) themselves, as long as these EOMs describe the movement of all of the state variables to be estimated. In practice, this requires defining a hyperdual space whose size is equal to the number of variables in the state to be estimated. For example, if estimating the Cartesian state a spacecraft and a maneuver magnitude at a reference point during an orbital determination arc, building a sevendimensional dual space will return the result of the EOMs and the partials matrix computed at that point. An open source example of the building of such hyperdual space has been implemented as a test case in \textit{dual_num}, a thorough dual number library in a programming language called Rust.

That's me.
\o/
↩