Keyboard shortcuts

Press ← or → to navigate between chapters

Press S or / to search in the book

Press ? to show this help

Press Esc to hide this help

Numerical validation

The evidence hierarchy for each numerical API is:

  1. direct unit cases and invalid-input tests;
  2. golden values generated by the pinned C++ implementation;
  3. independent mathematical properties or external reference data;
  4. Rust/Python cross-interface parity;
  5. equivalent release-mode benchmarks.

Golden files use schema-versioned JSON. Floating-point values are encoded as C99 hexadecimal strings, with NaN, +Infinity, and -Infinity used only when recording upstream non-finite behavior. This preserves every finite binary64 value without relying on a JSON parser’s decimal conversion.

The initial C++ baseline was built in release mode with GCC 14.3.0. All 29 upstream C++ test executables passed on 2026-07-25. Benchmark executables were built but their results are not presented as Rust comparisons until equivalent Rust workloads exist.

Randomized oracle cases use a named PCG32 seed recorded in each data file. Discovered failures are promoted to fixed regression cases.

Phase 3 adds 9 epoch cases and 134 anomaly cases from the pinned C++ implementation. The anomaly set includes every conversion direction, elliptic and hyperbolic boundary regimes, angles outside one revolution, and 64 deterministic PCG32 solver samples. The Rust layer additionally rejects invalid calendars and true anomalies outside the physical hyperbolic asymptote instead of propagating non-finite values.

Phase 4 adds 206 direct element/Cartesian cases and 16 analytic Jacobians from the pinned implementation. Jacobian files explicitly record row-major output-by-input order. Independent tests cover 2,000 deterministic elliptic and hyperbolic state round trips, both equinoctial pole conventions, finite-difference derivatives, and inverse-Jacobian identities. NumPy batches are compared row-for-row with the scalar Python API.

Phase 5 adds 28 propagation cases spanning zero and negative duration, circular, elliptic, hyperbolic, near-parabolic, and many-period trajectories. Each case records both Lagrange-coefficient and universal-variable states plus the Lagrangian and Reynolds 6 by 6 STMs. The source file SHA-256 is 63e2340aa8e4f65a4ae831b7e1c8eb3a28930e2e1ed336d9324e0dbf44469e57. Independent tests check energy and angular-momentum conservation, time reversal, central finite differences, and STM composition. Scalar propagation uses only fixed-size stack values; the core path performs zero heap allocations. Batch APIs allocate exactly their returned output storage after copying Python-owned input before releasing the GIL.

Phase 6 records transfer, encoding, flyby, MIMA, and 13 Lambert solutions across three geometries in phase6-v1.json (SHA-256 d2a81311264ecc94a6caee15854273abb4df5f32f682fe5340ffeb247e3dd3e3). Lambert ordering is zero revolution followed by left/right pairs. Every returned branch is independently propagated to its requested endpoint; encoding pairs round-trip, and the flyby Jacobian is checked by central differences. MIMA2 is checked against the published upstream reference case.

Phase 7 adds six Keplerian states over negative, reference, near-reference, and long-span MJD2000 epochs. phase7-v1.json has SHA-256 7eda8fb03796ab7d70cdb0e94c422773422ed682aff9310cd06512c0f6396a54. Independent tests cover period recurrence, all supported element representations, ordered batches, explicit unsupported capabilities, and concurrent read-only evaluation through shared ownership.

Phase 8 adds true-anomaly elements, Cartesian states, and physical metadata for all eight JPL low-precision bodies at five epochs spanning the open 1800–2050 validity interval. phase8-v1.json contains 40 cases and has SHA-256 0a3408893b5c04fdfddb452408057faa92d38a46775290c32eb6d2393683e3da. Independent tests cover case-insensitive lookup, the exact supported-name set, ordered batches, safe-radius validation, and both excluded interval boundaries. The provider retains the source table’s slightly negative fitted Earth inclination near the ends of the interval while keeping the general public classical-element converter’s canonical [0, π] inclination contract.

Phase 9 adds 54 VSOP2013 states covering all nine bodies at six epochs around the 1890–2000 fit interval and the J2000/MJD2000 half-day offset. Two additional cases verify coefficient selection at the default 1e-5 and coarse 0.5 thresholds. phase9-v1.json has SHA-256 9d5af01df18acb17fcd1b2356af6f2cc08f1cf75e10864270d0c30732f8b00d8. The oracle uses pykep 3.0.1 linked to heyoka 7.10.0. At the embedded high-precision floor of 1e-9, the maximum observed Rust/C++ difference is 0.185 m in position and 6.6e-8 m/s in velocity. Independent checks cover case-insensitive names, feature reporting, threshold errors, clone determinism, ordered Python batches, and a build/test run without the optional coefficient feature.

Phase 10 selects a pure-Rust DOP853 backend through a pykep-owned facade. Decision tests compare Kepler states and the 6 by 6 STM to the independent analytic propagator, verify a parameter-sensitivity column by central differences, bound 100-orbit energy drift, and exercise a CR3BP close approach with Jacobi drift and backward reversal checks. Separate tests cover rejected steps, step exhaustion, bit-for-bit repeated solves, dense interpolation, terminal root location, malformed grids, non-finite values, and physical singularities. The selected final-state path retains no internal trajectory and performs no per-step heap allocation for fixed-size model states.

The candidate harness and matching C++ benchmark use the same six-state Kepler initial condition, final time, and 1e-12 scalar tolerances. Nominal and variational timings, allocation limitations, and the dense-output maximum-step caveat are recorded in ADR 0004. Those tests validate integration machinery independently of the production models added in Phase 11.

Phase 11 adds five sampled states and the final 6 by 6 STM for each of Kepler, CR3BP, and BCP in phase11-v1.json. The BCP case uses a nonzero Sun mass and a nonzero initial epoch, while the CR3BP case reproduces the representative upstream trajectory. The file has SHA-256 1c2b67eb203da62baf921db9f811b0ad0f763cbaa03f0ca80d6319870b0da807. The oracle uses pykep 3.0.1 and heyoka 7.10.0 at a requested tolerance of 1e-16; Rust comparison settings and achieved tolerances are documented in dynamics.md. Separate tests evaluate the source equations directly, verify the triangular equilibrium, preserve the Jacobi constant, make zero-Sun BCP reduce to CR3BP, check every state/parameter Jacobian column by central differences, and distinguish body singularities from solver failures.

Phase 12 records the final state and all first-order state/control variations for one upstream regression case from each of the four ZOH systems. phase12-v1.json has SHA-256 294bab93355628cd22991b83563c6f93deefc19dd93911dfc4b654a101764f22. Single-segment nominal tolerances range from 2e-12 to 3e-10; variation tolerances range from 2e-7 to 2e-5, reflecting the fixed-size numerical Jacobians documented in zero-order-hold.md. Independent tests cover exact switch ownership, malformed grids, zero-control reductions, manual segment-by-segment equivalence, forward/backward reversal, and activation of per-segment sensitivity columns.

Phase 13 records mass- and time-optimal Cartesian and modified-equinoctial states and all first-order variations with respect to the seven initial costates and upstream lambda0 variational argument. It also records the upstream dimensional 100-day equinoctial case. phase13-v1.json has SHA-256 a1276c4c35c7ad60b481d81d69f8df64b542bad1d5cc6571738db820cb0e2c3d. Normalized nominal trajectories use a scaled 3e-10 bound, the dimensional case uses 3e-9, and variations use 2e-4 to account for the generic centered numerical Jacobians. Independent tests cover minimized-Hamiltonian conservation, propagated central differences, output-by-input Jacobian orientation, explicit zero-primer errors, control normalization, and canonical Hamiltonian agreement across the analytic coordinate/costate transform.

Phase 14 records fixed-duration Sims–Flanagan mismatch and throttle constraints, all three analytic mismatch-Jacobian groups, and throttle Jacobians for one physical five-segment case and normalized four-segment cases at cuts zero, one half, and one. It also records equal and irregular direct durations for the alpha variant. phase14-v1.json has SHA-256 0bd6adddc72d850f1de5e4ab95a436128d9b8181f8de2b3bbbc67300e004d542. The Rust analytic gradients agree with the pinned C++ values to scaled bounds of 2e-10 to 3e-10 and with independent scale-adjusted central differences. Separate tests cover one segment, odd/even splits, zero and unit-limit throttle, all cut boundaries, normalized weights, invalid propulsion/gravity/mass values, dimension mismatches, and non-finite inputs.

Phase 15 records generic ZOH-leg mismatches and all four Jacobian groups for Kepler, CR3BP, modified-equinoctial, and ideal solar-sail dynamics. phase15-v1.json has SHA-256 9d8d424c2af10ceee497f74f62acb30c53924f78f7be8a528dbe44bf14767935. Nominal mismatches agree with the pinned C++/heyoka oracle within a scaled 3e-9 bound; endpoint, chronological-control, and time-grid derivatives agree within 3e-5, reflecting the fixed-size numerical model Jacobians. Every Kepler derivative column is also checked against independent scale-adjusted central differences. Separate tests cover cut zero and one, strict time grids, dimensions, finite values, solver options, state histories, ordered batches, and maximum-step exhaustion with direction, segment index, and time interval in the reported failure.

Phase 16 audits the complete Python surface against the pinned upstream exports and records every equivalent, renamed, intentionally different, deferred, and unsupported area in python-migration.md. Pytest checks every runtime export against the shipped stub, including parameter names/order, defaults, return-annotation presence, and class-member documentation, then exercises strided/read-only arrays, wrong dtype and rank, wrong shape, NaN/infinity, typed exceptions, repeatability, copied constructor inputs, and concurrent reuse. Mypy checks representative code against the packaged stub. A clean wheel matrix installs and imports CPython 3.11–3.13 wheels on Linux, macOS, and Windows; local shared-library inspection verifies that the extension has no C++ runtime dependency.

Phase 17 adds deterministic Rust/Python pairs for epochs and anomalies, elements and propagation, Lambert branches, ephemerides, gravity assists, Sims–Flanagan gradients, and CR3BP/ZOH dynamics, plus a NumPy batch example. Every Rust binary is compiled under workspace test and clippy gates and run locally in release mode. Pytest discovers and runs every Python script against the installed release extension. Ten Rustdoc examples cover every major module landing page and representative Lambert/ZOH types and compile with warnings denied.

Phase 18 adds a protocol-matched 100-sample Rust/C++ performance distribution, coarse CI regression thresholds, batch scaling, and allocation/cache/ vectorization profiling before any optimization. Miri runs suitable structural and parser tests, Valgrind checks the release harness, and bounded libFuzzer/AddressSanitizer campaigns cover epoch parsing, element conversion, Lambert inputs, and Reynolds-STM overflow boundaries. stabilization.md records exact results, limits, tool constraints, and the absence of algorithm changes.