astrokit

A highly-performant astrodynamics engine, written in Rust. Force models, propagation, reference frames, spacecraft geometry and measurement models — built on ordinary-diffeq for integration, and generic over the scalar type so that any propagation you can solve you can also differentiate.

Status

The crate is still being built to a milestone plan:

• M1 — the core subset: two-body and spherical-harmonic gravity, drag, third-body, cannonball SRP, GCRF/ITRF frames, ground stations and the four measurement types, with a validation suite against Orekit reference trajectories.
• M2 — general-purpose completion: higher-degree gravity fields, more density models, box-wing SRP and eclipses, maneuvers, covariance propagation, relative motion, TLE and CCSDS I/O.

The tables under Validation and Performance are the acceptance criteria.

Conventions

• Everything is SI
• Time is typed, using hifitime
• US English spelling

Foundations

These are depended on, not reimplemented.

Crate	License	Provides
hifitime 4.3	MPL-2.0	TAI/TT/ET/TDB/UTC/GPS time scales, leap seconds
anise 0.10	MPL-2.0	Reference frames, planetary ephemerides, SPICE-format kernels
ordinary-diffeq 0.8	MIT	Adaptive integration, dense output, callbacks, autodiff

Validation

Before we can make a highly-performant propagator, we have to make a correct propagator.

1. Orekit reference fixtures. A committed, rerunnable generation script produces reference trajectories; the data lands under tests/fixtures/ with the Orekit version, force-model configuration and generation date recorded beside it, and Rust tests assert agreement within declared tolerances. Observations are still to come.
2. Analytic invariants.: ie, "does a two-body orbit close back in on itself after one-period?" and similar conservation-tests, meant to run fast as unit tests to ensure that nothing has broken.
3. ANISE cross-checks. Frame transformations and ephemeris lookups against ANISE's own validated values, confirming that we are calling it correctly.

Every fixture comparison matches gravity degree and order, density model, third-body ephemeris source and integrator tolerances; anything that cannot be matched is documented and its contribution bounded.

Against Orekit

One force at a time over a one-day low Earth arc, worst disagreement across all 1441 states, against Orekit 13.1.8. Each force is isolated so that a disagreement is attributable; the figures are read against the point-mass row, where there is no model left to mismatch and what is measured is the comparison itself — frames, time scale, units and mu.

Case	Target	Achieved	Test
Point mass	1 m, 1 mm/s	0.77 mm, 8.2e-7 m/s	point_mass_arc
Sun and Moon, over a point mass	1 m, 1 mm/s	0.84 mm, 8.5e-7 m/s	sun_and_moon
Cannonball radiation pressure	1 m, 1 mm/s	1.00 mm, 1.0e-6 m/s	radiation_pressure
EGM96 8x8	1 m, 1 mm/s	1.09 cm, 9.1e-6 m/s	low_degree_field
EGM96 20x20	1 m, 1 mm/s	1.10 cm, 9.0e-6 m/s	higher_degree_field
Drag, NRLMSISE-00	1e-3% of the drag effect	2.7 cm, 3.1e-5 m/s — 0.0001%	msis_drag
6-day LEO arc	to be set from the 1-day result	not yet measured	—
GEO arc, including a time-varying-area change	to be set	not yet measured	—
RA/DEC, range, range-rate, az/el against Orekit	to be set	not yet measured	—

Tests are in tests/orekit_comparison.rs; run them with --nocapture to see the figures above printed.

The drag row is bounded as a fraction of what drag does over the arc rather than in meters, because the meters depend on the arc and the fraction does not. Configured identically, two independent NRLMSISE-00 implementations agree to four parts per million of the drag effect.

Harris-Priester is deliberately not compared against either tool. This crate departs from the published interpolation on purpose — it fits a monotone cubic to the logarithm of the density so the altitude derivative exists, where the others interpolate the table exponentially — so a comparison measures that choice rather than correctness. The model itself is tested against the published form in src/atmosphere/harris_priester.rs.

Against nyx-space

A second, independent reference. nyx reads the same ANISE kernels this crate does, so where the Orekit comparison is limited by two realizations of the inertial frame being 23 mas apart, here both sides use literally the same frame — and for gravity the same ITRF93 out of the same binary PCK. That removes the term that dominates the Orekit gravity rows, and the agreement is two to three orders of magnitude tighter.

Sampling differs too. Orekit's fixtures are a 60-second grid; nyx's are its own accepted integrator steps, because its trajectory interpolation is a five-point Hermite fitted after the fact rather than a continuous extension of the method, and that interpolation is its accuracy floor. On the point-mass arc its stored states agree with Kepler to 2.7e-6 m while the same trajectory read on a 60-second grid is 8.1e-4 m.

Case	Target	Achieved	Test
Point mass	1 m, 1 mm/s	0.218 mm, 2.4e-7 m/s	point_mass_arc
Sun and Moon, over a point mass	1 m, 1 mm/s	0.197 mm, 2.1e-7 m/s	sun_and_moon
EGM96 8x8	1 m, 1 mm/s	0.102 mm, 1.1e-7 m/s	low_degree_field
EGM96 20x20	1 m, 1 mm/s	0.021 mm, 2.2e-8 m/s	higher_degree_field
Cannonball radiation pressure	1 m, 1 mm/s	0.491 mm, 5.3e-7 m/s	radiation_pressure
Drag, NRLMSISE-00	1% of the drag effect	8.55 m, 9.7e-3 m/s — 0.029%	msis_drag

Tests are in tests/nyx_comparison.rs, fixtures from scripts/generate-nyx-fixtures.sh. nyx is AGPL-3.0-or-later, so it is run offline to produce data and is not a dependency of anything here.

The drag row is the only one where the two sides disagree by more than a millimeter, and the reason is known. Both tools are given numerically identical space weather — the fixture's drivers are written out of the same record this crate reads, and nyx takes the 81-day centered average from a file column, so it carries this crate's own value rather than a separately published one. With the inputs identical the densities still differ by 2.6e-4. Run against case 0 of the seventeen standard NRLMSISE-00 cases, which pin the NRL reference to sixteen digits, this crate reproduces the published total mass density exactly and nyx is 1.08e-3 from it. So that 0.029% is nyx departing from the reference.

One element is unmatched and it is this crate's side: nyx measures geodetic altitude against the planetary constants kernel's Earth and this crate against WGS84, putting the same spacecraft 0.400 m apart. That is worth about 7e-6 of the density, two orders below the model difference, and it cannot be removed from the nyx side — ANISE replaces a frame's ellipsoid from the loaded kernel.

Analytic invariants

Case	Target	Achieved	Test
Two-body closure over one period	integrator tolerance	9.8e-7 m, 1.1e-9 m/s	orbit_closes_over_one_period
Two-body energy and angular momentum, five periods	integrator tolerance	1.1e-13, 5.3e-14 relative	conserves_energy_and_angular_momentum
J2 secular nodal regression and apsidal rotation	1%	3.8e-3, 6.2e-3 relative	j2_secular_rates_match_the_closed_form

The J2 residual is systematic rather than noise: the closed form is first order and describes mean elements where the fit is on osculating ones, and the field's symmetry axis is the true pole while the inclination is measured in J2000. It is a check against a sign, normalization or coefficient error, all of which would miss by far more than a percent. The tight comparison is the Orekit one above.

Performance

Measured against nyx-space and Orekit on the same benchmarks with matched force models, on the same hardware. Neither publishes usable figures, so the references are measured here rather than cited.

The scenarios are the ones the fixtures above were generated from, so the trajectories being timed are known to agree. One day of arc, steps and milliseconds, measured 2026-10-02. Orekit runs under a JVM warmed ten times before anything is recorded, and the minimum of thirty runs is reported rather than the mean.

Benchmark	astrokit	nyx-space	Orekit
Two-body	1752 steps, 0.49 ms	2151 steps, 7.06 ms	1932 steps, 4.60 ms
EGM96 8x8	1968 steps, 31.9 ms	2158 steps, 147.5 ms	1938 steps, 34.2 ms
EGM96 20x20	2473 steps, 54.5 ms	2158 steps, 193.6 ms	1940 steps, 52.9 ms
Sun and Moon	1794 steps, 35.9 ms	2152 steps, 79.7 ms	1932 steps, 36.1 ms
Radiation pressure	2067 steps, 28.0 ms	2570 steps, 166.9 ms	2035 steps, 133.6 ms
Drag, NRLMSISE-00	1859 steps, 52.0 ms	2452 steps, 385.0 ms	2092 steps, 357.2 ms
Ensemble of N states, scaling vs. cores	not yet measured	—	—

Step counts are reported alongside the times because the tolerances are not literally equal: Orekit is asked for a nanometer of position accuracy and this runs at Accuracy::High. The three do comparable work — within about a third of each other on every row — so a ratio of times is close to a ratio of per-step cost. A future scenario whose step counts diverge needs them brought back together before its ratio means anything.

The radiation row is not like for like however the steps come out. Orekit's radiation pressure registers eclipse event detectors, so its integrator locates the shadow boundaries by root finding, and this crate integrates straight through them. That is work Orekit does and this does not.

nyx is timed the same way, by the minimum of thirty runs after warmup, and against a control that skips building the trajectory so the harness can be told apart from the propagation. The control is within noise on every row, so these are the propagations.

One structural difference is worth knowing before reading the two-body row. nyx's integrated state is Const<90> — nine components plus an eighty-one element state transition matrix — and it is carried on every step whether or not an STM was asked for, where this crate integrates the six a two-body problem has. On a trivial right-hand side that vector arithmetic is most of the work, which is why the gap is widest where the physics is cheapest and narrows as the force model gets expensive.

Reproduce with cargo bench --bench reference_scenarios, scripts/generate-orekit-fixtures.sh --time and scripts/generate-nyx-fixtures.sh --time.

Kernels

Ephemerides and Earth orientation come from NAIF kernels, which are too large to commit. data/kernels.tsv pins each one by URL and SHA-256:

scripts/fetch-kernels.sh            # into data/kernels
scripts/fetch-kernels.sh /somewhere # or wherever you keep them

The script verifies checksums and leaves an already-correct file alone. Point ASTROKIT_KERNELS at the directory if it is not data/kernels. Tests that need kernels say so and skip when they are absent.

earth_latest_high_prec.bpc is deliberately not used: a result that cites "the latest Earth orientation kernel" cannot be reproduced, so the dated snapshot covering the same span is pinned instead.

Examples

cargo run --example two_body   # one orbital period of a point mass
cargo run --example attitude   # orbit and attitude together, on a tumbling box
cargo run --example full       # a day of EGM96 gravity, with the ground track

full needs kernels. It takes about 6 seconds unoptimized and 3 with --release, because it propagates the day seven times: once for the trajectory, once for a bare point mass, and once per force with that force left out, which is how it prices each of them.

Development

Enable the pre-commit hook once per clone:

git config core.hooksPath .githooks

It runs cargo fmt --check, cargo clippy --all-targets --all-features -- -D warnings, and a US English spelling check. git commit --no-verify skips it. Tests and docs are left out to keep it under a few seconds; run cargo test yourself.

License

MIT. See LICENSE.