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