Skip to main content

astrodyn_quantities/
integ_origin.rs

1//! Integration-frame origin and the typed shift to root-inertial.
2//!
3//! `Position<IntegrationFrame>` lives in a body's integration frame — a
4//! non-rotating frame whose origin generally differs from the root
5//! inertial origin (e.g. integrating a vehicle in `Earth.inertial` when
6//! the root frame is `SSB.inertial`).
7//!
8//! Apply this shift only when the consumer requires root-inertial
9//! coordinates — i.e. when it mixes body state with root-inertial source
10//! positions (Sun, Moon, gravity sources). Those *shift sites* are:
11//! gravity, relativistic corrections, SRP, solar beta, earth lighting.
12//!
13//! Do **not** shift for atmosphere, drag velocity, LVLH, geodetic, or
14//! orbital elements: those operate within a single planet's inertial
15//! frame, which is *the body's integration frame* in realistic configs,
16//! so `body.trans.position` is already the correct input. They take
17//! [`Position<PlanetInertial<P>>`](crate::frame::PlanetInertial) at the
18//! typed-sibling API.
19//!
20//! See issue #255 / `RF.10` in `docs/JEOD_invariants.md` for the
21//! shift-vs-no-shift split that motivated this module.
22
23use crate::aliases::{Position, Velocity};
24use crate::frame::{IntegrationFrame, RootInertial};
25
26/// The integration-frame origin: position and velocity of the
27/// integration-frame origin expressed in root-inertial coordinates.
28///
29/// Construct from the runner's frame-tree state. `IntegOrigin::zero()` is
30/// used when the body integrates in the root frame, so the shift is a
31/// no-op.
32///
33/// # RF.10 structural guard
34///
35/// `Position<IntegrationFrame>` and `Position<RootInertial>` are
36/// kind-distinct frames: the only Sub/Add impls require
37/// `CompatibleFrames<F, F>`, so subtracting one from the other is a
38/// compile error. The reviewer-flagged bug shape from PR #258 — mixing
39/// integration-frame body state with root-inertial source positions —
40/// cannot recur silently.
41///
42/// **`Position<IntegrationFrame> - Position<RootInertial>` does not compile:**
43///
44/// ```compile_fail
45/// use astrodyn_quantities::prelude::*;
46/// let body: Position<IntegrationFrame> = Position::zero();
47/// let sun: Position<RootInertial> = Position::zero();
48/// // Frames mismatch — `CompatibleFrames<IntegrationFrame, RootInertial>` unimplemented:
49/// let _bug = body - sun;
50/// ```
51///
52/// **Adding them does not compile either:**
53///
54/// ```compile_fail
55/// use astrodyn_quantities::prelude::*;
56/// let body: Position<IntegrationFrame> = Position::zero();
57/// let sun: Position<RootInertial> = Position::zero();
58/// let _bug = body + sun;
59/// ```
60///
61/// **Cross-assignment refuses:**
62///
63/// ```compile_fail
64/// use astrodyn_quantities::prelude::*;
65/// let body: Position<IntegrationFrame> = Position::zero();
66/// let _bug: Position<RootInertial> = body;   // not the same type
67/// ```
68///
69/// The only way through is the typed shift —
70/// [`shift_position`](IntegOrigin::shift_position) for step-constant or
71/// [`shift_position_at_stage`](IntegOrigin::shift_position_at_stage)
72/// for stage-time-interpolated:
73///
74/// ```
75/// use astrodyn_quantities::prelude::*;
76/// let o = IntegOrigin::zero();
77/// let body: Position<IntegrationFrame> = Position::zero();
78/// let body_root: Position<RootInertial> = o.shift_position(body);
79/// let sun: Position<RootInertial> = Position::zero();
80/// let sun_to_vehicle: Position<RootInertial> = body_root - sun;   // OK
81/// # let _ = sun_to_vehicle;
82/// ```
83#[derive(Debug, Clone, Copy, Default)]
84pub struct IntegOrigin {
85    /// Position of the integration-frame origin in root-inertial coordinates.
86    pub position: Position<RootInertial>,
87    /// Velocity of the integration-frame origin in root-inertial coordinates.
88    pub velocity: Velocity<RootInertial>,
89}
90
91impl IntegOrigin {
92    /// The zero offset — used when the integration frame coincides with
93    /// the root inertial frame.
94    #[inline]
95    pub fn zero() -> Self {
96        Self {
97            position: Position::<RootInertial>::zero(),
98            velocity: Velocity::<RootInertial>::zero(),
99        }
100    }
101
102    /// Shift a position from the integration frame to root inertial by
103    /// adding the origin's root-inertial position.
104    #[inline]
105    pub fn shift_position(&self, p: Position<IntegrationFrame>) -> Position<RootInertial> {
106        Position::<RootInertial>::from_raw_si(p.raw_si() + self.position.raw_si())
107    }
108
109    /// Shift a velocity from the integration frame to root inertial by
110    /// adding the origin's root-inertial velocity.
111    #[inline]
112    pub fn shift_velocity(&self, v: Velocity<IntegrationFrame>) -> Velocity<RootInertial> {
113        Velocity::<RootInertial>::from_raw_si(v.raw_si() + self.velocity.raw_si())
114    }
115
116    /// Stage-time-interpolated shift to root inertial.
117    ///
118    /// Used inside RK4 derivative closures where the intermediate
119    /// position is sampled at `time_frac * dt` into the step. The
120    /// integration-frame origin itself moves under its own velocity, so
121    /// the stage-time shift is `origin.position + origin.velocity *
122    /// stage_dt` (linear interpolation; matches the arithmetic used by
123    /// the runner's gravity stage closure, RF.10). For
124    /// integration-frames at rest in root, this collapses to the
125    /// step-constant `shift_position`.
126    #[inline]
127    pub fn shift_position_at_stage(
128        &self,
129        p: Position<IntegrationFrame>,
130        stage_dt: f64,
131    ) -> Position<RootInertial> {
132        let stage_origin = self.position.raw_si() + self.velocity.raw_si() * stage_dt;
133        Position::<RootInertial>::from_raw_si(p.raw_si() + stage_origin)
134    }
135}
136
137#[cfg(test)]
138mod tests {
139    use super::*;
140    use glam::DVec3;
141
142    #[test]
143    fn zero_origin_is_identity_for_position() {
144        let p = Position::<IntegrationFrame>::from_raw_si(DVec3::new(7e6, 0.0, 0.0));
145        let shifted = IntegOrigin::zero().shift_position(p);
146        assert_eq!(shifted.raw_si(), p.raw_si());
147    }
148
149    #[test]
150    fn zero_origin_is_identity_for_velocity() {
151        let v = Velocity::<IntegrationFrame>::from_raw_si(DVec3::new(0.0, 7500.0, 0.0));
152        let shifted = IntegOrigin::zero().shift_velocity(v);
153        assert_eq!(shifted.raw_si(), v.raw_si());
154    }
155
156    #[test]
157    fn nonzero_origin_adds_to_position() {
158        let o = IntegOrigin {
159            position: Position::<RootInertial>::from_raw_si(DVec3::new(1.5e11, 0.0, 0.0)),
160            velocity: Velocity::<RootInertial>::zero(),
161        };
162        let p_integ = Position::<IntegrationFrame>::from_raw_si(DVec3::new(7e6, 0.0, 0.0));
163        let p_root = o.shift_position(p_integ);
164        assert_eq!(p_root.raw_si(), DVec3::new(1.5e11 + 7e6, 0.0, 0.0));
165    }
166
167    #[test]
168    fn nonzero_origin_adds_to_velocity() {
169        let o = IntegOrigin {
170            position: Position::<RootInertial>::zero(),
171            velocity: Velocity::<RootInertial>::from_raw_si(DVec3::new(0.0, 30_000.0, 0.0)),
172        };
173        let v_integ = Velocity::<IntegrationFrame>::from_raw_si(DVec3::new(0.0, 7500.0, 0.0));
174        let v_root = o.shift_velocity(v_integ);
175        assert_eq!(v_root.raw_si(), DVec3::new(0.0, 30_000.0 + 7500.0, 0.0));
176    }
177}