Skip to main content

sidereon_core/astro/
mod.rs

1//! Numerical astrodynamics engine for orbit propagation, force models, and
2//! future flight-dynamics primitives.
3//!
4//! Current scope:
5//!
6//! - inertial Cartesian state representation
7//! - two-body gravity and J2 perturbation
8//! - fixed-step RK4
9//! - adaptive Dormand-Prince 5(4) (`DP54`)
10//! - propagation results with step statistics
11//!
12//! Planned future work includes dense output, event handling, richer propagation
13//! contexts, additional force models, covariance propagation, estimation, and
14//! maneuver support.
15
16pub mod almanac;
17pub mod angles;
18pub mod anomaly;
19pub mod apparent;
20pub mod atmosphere;
21pub mod bodies;
22pub mod cdm;
23pub mod conjunction;
24pub mod constants;
25pub mod covariance;
26pub mod coverage;
27pub mod data;
28pub mod doppler;
29pub mod elements;
30pub mod equinoctial;
31pub mod error;
32pub mod events;
33pub mod forces;
34pub mod frames;
35pub mod integrators;
36pub mod iod;
37pub mod lambert;
38pub mod math;
39pub mod ndm;
40pub mod observation;
41pub mod oem;
42pub mod omm;
43pub mod opm;
44pub mod passes;
45pub mod propagator;
46pub mod relative;
47pub mod rf;
48pub mod sgp4;
49pub mod spk;
50pub mod state;
51pub mod tca;
52pub mod time;
53pub mod tle;
54pub mod tolerances;
55pub mod xml;
56
57pub use spk::{
58    DafByteOrder, DafFileRecord, DafSpk, Spk, SpkError, SpkSegmentDescriptor, SpkState,
59    SpkStateVector,
60};
61
62#[cfg(all(feature = "sgp4-debug-oracle", sgp4_oracle_built))]
63#[doc(hidden)]
64pub mod sgp4_cpp_oracle {
65    //! Test-only oracle bridge to the Vallado C++ implementation.
66    //! Compiled in only when the `sgp4-debug-oracle` feature is on and the
67    //! development-only C++ oracle sources were found by the build script.
68    //! Not part of the public API.
69
70    use std::os::raw::{c_char, c_double, c_int};
71
72    pub const CPP_DUMP_DOUBLE_COUNT: usize = 112;
73    pub const CPP_DUMP_INT_COUNT: usize = 5;
74
75    extern "C" {
76        pub fn cpp_sgp4init_dump(
77            satnum: *const c_char,
78            opsmode: c_char,
79            epoch_sgp4: c_double,
80            bstar: c_double,
81            ndot: c_double,
82            nddot: c_double,
83            ecco: c_double,
84            argpo: c_double,
85            inclo: c_double,
86            mo: c_double,
87            no_kozai: c_double,
88            nodeo: c_double,
89            epochyr: c_int,
90            epochdays: c_double,
91            jdsatepoch: c_double,
92            jdsatepoch_frac: c_double,
93            double_out: *mut c_double,
94            int_out: *mut c_int,
95        ) -> c_int;
96
97        pub fn cpp_sgp4_step(
98            satnum: *const c_char,
99            opsmode: c_char,
100            epoch_sgp4: c_double,
101            bstar: c_double,
102            ndot: c_double,
103            nddot: c_double,
104            ecco: c_double,
105            argpo: c_double,
106            inclo: c_double,
107            mo: c_double,
108            no_kozai: c_double,
109            nodeo: c_double,
110            epochyr: c_int,
111            epochdays: c_double,
112            jdsatepoch: c_double,
113            jdsatepoch_frac: c_double,
114            tsince: c_double,
115            r_out: *mut c_double,
116            v_out: *mut c_double,
117        ) -> c_int;
118    }
119
120    /// Force-reference the C symbols so the linker pulls in the static lib.
121    /// Without this, the rlib has no use of the symbols and the linker
122    /// strips the entire archive when compiling integration tests.
123    #[doc(hidden)]
124    pub fn force_link_oracle() -> usize {
125        let init_dump = cpp_sgp4init_dump as *const ();
126        let step = cpp_sgp4_step as *const ();
127
128        init_dump as usize ^ step as usize
129    }
130}
131
132#[cfg(all(feature = "sgp4-debug-oracle", sgp4_oracle_built))]
133pub use sgp4_cpp_oracle::cpp_sgp4_step;
134
135pub use anomaly::{
136    eccentric_to_mean, eccentric_to_true, mean_to_eccentric, mean_to_true, propagate_kepler,
137    solve_kepler, true_to_eccentric, true_to_mean, AnomalyError, KeplerSolution,
138};
139pub use elements::{coe2rv, rv2coe, ClassicalElements, ElementsError, OrbitType};
140pub use equinoctial::{
141    coe2eq, coe2mee, eq2coe, eq2mee, eq2rv, mee2coe, mee2eq, mee2rv, rv2eq, rv2mee,
142    EquinoctialElements, EquinoctialError, ModifiedEquinoctialElements, RetrogradeFactor,
143};
144pub use error::PropagationError;
145pub use state::CartesianState;
146pub use time::Time;
147
148#[cfg(test)]
149mod tests {
150    use super::*;
151    use crate::astro::forces::TwoBodyGravity;
152    use crate::astro::integrators::{Integrator, DP54};
153    use crate::astro::propagator::{api::IntegratorOptions, OrbitalDynamics, PropagationContext};
154    use nalgebra::Vector3;
155
156    #[test]
157    fn test_two_body_dp54_precision() {
158        let r_mag: f64 = 7000.0;
159        let mu: f64 = 398600.4418;
160        let v_mag: f64 = (mu / r_mag).sqrt();
161        let initial_state = CartesianState {
162            epoch_tdb_seconds: 0.0,
163            position_km: Vector3::new(r_mag, 0.0, 0.0),
164            velocity_km_s: Vector3::new(0.0, v_mag, 0.0),
165        };
166
167        let force = TwoBodyGravity::default();
168        let dynamics = OrbitalDynamics {
169            force_model: &force,
170        };
171        let integrator = DP54;
172        let ctx = PropagationContext::default();
173        let opts = IntegratorOptions {
174            abs_tol: 1e-12,
175            rel_tol: 1e-12,
176            initial_step: 1.0,
177            min_step: 1e-15,
178            ..IntegratorOptions::default()
179        };
180
181        let period = 2.0 * std::f64::consts::PI * (r_mag.powi(3) / mu).sqrt();
182        let result = integrator
183            .propagate(initial_state, period, &dynamics, &ctx, &opts)
184            .unwrap();
185
186        let final_pos = result.final_state.position_km;
187        let final_vel = result.final_state.velocity_km_s;
188
189        // Oracle 1: Return to start precision (Sub-millimeter)
190        assert!(
191            (final_pos.x - r_mag).abs() < 1e-7,
192            "Position X error too large: {}",
193            (final_pos.x - r_mag).abs()
194        );
195        assert!(
196            final_pos.y.abs() < 1e-7,
197            "Position Y error too large: {}",
198            final_pos.y.abs()
199        );
200
201        // Oracle 2: Energy conservation (Specific mechanical energy)
202        let initial_energy = v_mag.powi(2) / 2.0 - mu / r_mag;
203        let final_v_mag = final_vel.norm();
204        let final_r_mag = final_pos.norm();
205        let final_energy = final_v_mag.powi(2) / 2.0 - mu / final_r_mag;
206        assert!(
207            (final_energy - initial_energy).abs() < 1e-10,
208            "Energy conservation failure: {}",
209            (final_energy - initial_energy).abs()
210        );
211    }
212}