sidereon_core/astro/forces/
j2.rs

1use crate::astro::constants::{J2_EARTH, MU_EARTH, RE_EARTH};
2use crate::astro::error::PropagationError;
3use crate::astro::forces::r#trait::ForceModel;
4use crate::astro::propagator::api::PropagationContext;
5use crate::astro::state::CartesianState;
6use nalgebra::Vector3;
7
8pub struct J2Gravity {
9    pub mu: f64,
10    pub re: f64,
11    pub j2: f64,
12}
13
14impl Default for J2Gravity {
15    fn default() -> Self {
16        Self {
17            mu: MU_EARTH,
18            re: RE_EARTH,
19            j2: J2_EARTH,
20        }
21    }
22}
23
24impl ForceModel for J2Gravity {
25    fn acceleration(
26        &self,
27        state: &CartesianState,
28        _ctx: &PropagationContext,
29    ) -> Result<Vector3<f64>, PropagationError> {
30        let r_mag2 = state.position_km.norm_squared();
31        if r_mag2 == 0.0 {
32            return Err(PropagationError::NumericalFailure(
33                "Zero position magnitude".to_string(),
34            ));
35        }
36        let r_mag = r_mag2.sqrt();
37
38        let f = 1.5 * self.j2 * (self.mu / r_mag2) * (self.re / r_mag).powi(2);
39        let z_r2 = (state.position_km.z * state.position_km.z) / r_mag2;
40
41        Ok(Vector3::new(
42            f * (state.position_km.x / r_mag) * (5.0 * z_r2 - 1.0),
43            f * (state.position_km.y / r_mag) * (5.0 * z_r2 - 1.0),
44            f * (state.position_km.z / r_mag) * (5.0 * z_r2 - 3.0),
45        ))
46    }
47}
48
49#[cfg(test)]
50mod tests {
51    use super::*;
52    use crate::astro::propagator::api::PropagationContext;
53    use crate::astro::state::CartesianState;
54
55    #[test]
56    fn acceleration_matches_orbis_force_wrapper_bits() {
57        let state = CartesianState::new(0.0, [7000.0, -1210.0, 1300.0], [0.0, 0.0, 0.0]);
58        let acceleration = J2Gravity::default()
59            .acceleration(&state, &PropagationContext::default())
60            .unwrap();
61
62        assert_eq!(acceleration.x.to_bits(), 13_754_131_348_549_160_135);
63        assert_eq!(acceleration.y.to_bits(), 4_519_025_615_523_880_849);
64        assert_eq!(acceleration.z.to_bits(), 13_750_824_904_549_515_386);
65    }
66}
sidereon_core/astro/forces/j2.rs

sidereon_core/astro/forces/
j2.rs
